Topic: Can't link openblas or lapack

I can't link to blas or lapack. In fact, I can't even find lapack or blas.

2 (edited by guitarfdm Today 13:45:01)

Re: Can't link openblas or lapack

I tried a reinstall as admin with virus software turned off, and it is still not here. There is an openblas_config.h in the C:\Program Files (x86)\Simply Fortran 3\mingw-w64\x86_64-w64-mingw32\include but no
libopenblas.a, libopenblas.dll.a, openblas.dll.

Microsoft Windows 11 Professional (x64) Build 26200.8655 (25H2)

It tries to auto-link and fails.

here is the compiler output:
Generating fortran-library.dll
build\burgMemcof.o: In function `linRegMA':
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
C:\Users\Admin\Fortran\Burg/./burgMemcof.f90:84: undefined reference to `dgels'
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
C:\Users\Admin\Fortran\Burg/./burgMemcof.f90:91: undefined reference to `dgels'
collect2.exe: error: ld returned 1 exit status
Error: Last command making (fortran-library.dll) returned a bad status
Error: Make execution terminatedcompiler output.

Re: Can't link openblas or lapack

Here is the module:

module burgMemcof
    use iso_c_binding, only: C_DOUBLE, C_INT
    implicit none

   
!    interface
!        subroutine dgels(trans, m, n, nrhs, a, lda, b, ldb, work, lwork, info) bind(C, name="dgels_")
!            ! Pull the types directly into the interface block itself
!            use, intrinsic :: iso_c_binding, only: C_CHAR, C_INT, C_DOUBLE
!            implicit none
!           
!            character(kind=C_CHAR) :: trans
!            integer(kind=C_INT)    :: m, n, nrhs, lda, ldb, lwork, info
!            real(kind=C_DOUBLE)    :: a(*), b(*), work(*)
!        end subroutine dgels
!    end interface

   
contains

subroutine linRegMA(residual, arCoefNum, maCoefNum, sampleSize, linRegEst) BIND(C, NAME="linRegMA")
    use iso_c_binding, only: C_DOUBLE, C_INT
    implicit none
    !GCC$ ATTRIBUTES DLLEXPORT :: linRegMA

    ! Inputs
    real(C_DOUBLE), dimension(*), intent(in) :: residual
    integer(C_INT), intent(in) :: arCoefNum, maCoefNum, sampleSize

    ! Output
    real(C_DOUBLE), dimension(*), intent(out) :: linRegEst

    ! Local variables
    integer(C_INT) :: m, n
    integer(C_INT) :: i, j
    real(C_DOUBLE), allocatable :: designMatrix(:,:), targetSet(:)
    real(C_DOUBLE), allocatable :: work(:)
    integer(C_INT) :: lwork, info
   
    external dgels
    ! -----------------------------
    ! Compute dimensions
    ! -----------------------------
    m = sampleSize - (arCoefNum + maCoefNum)
    n = maCoefNum + 1      ! +1 for intercept column

    ! -----------------------------
    ! Allocate arrays
    ! -----------------------------
    allocate(designMatrix(m, n))
    allocate(targetSet(m))

    ! -----------------------------
    ! Build targetSet
    ! Mathematica: residual[[arCoefNum + maCoefNum + 1 ;; sampleSize]]
    ! -----------------------------
    do i = 1, m
        targetSet(i) = residual(arCoefNum + maCoefNum + i)
    end do

    ! -----------------------------
    ! Build design matrix
    ! Column 1 = intercept = 1.0
    ! Columns 2..n = lagged residual  values
    ! -----------------------------
    designMatrix(:,1) = 1.0_C_DOUBLE

    do j = 1, maCoefNum
        do i = 1, m
            designMatrix(i, j+1) = residual(arCoefNum + j + (i-1))
        end do
    end do

    ! -----------------------------
    ! Solve Least Squares using DGELS
    ! A = designMatrix (m x n)
    ! B = targetSet (m x 1)
    ! DGELS overwrites B with solution coefficients
    ! -----------------------------
    lwork = -1
    allocate(work(1))

    ! Workspace query
    call dgels('N', m, n, 1, designMatrix, m, targetSet, m, work, lwork, info)

    lwork = int(work(1))
    deallocate(work)
    allocate(work(lwork))

    ! Actual solve
    call dgels('N', m, n, 1, designMatrix, m, targetSet, m, work, lwork, info)

    if (info /= 0) then
        print *, "DGELS failed with INFO=", info
    end if

    ! -----------------------------
    ! Extract coefficients
    ! DGELS stores solution in first n entries of targetSet
    ! -----------------------------
    linRegEst(1:n) = targetSet(1:n)

    ! Cleanup
    deallocate(designMatrix, targetSet, work)
end subroutine linRegMA


subroutine memcof(tseries, n, m, xms, d) BIND(C, NAME="memcof")
    implicit none
    !GCC$ ATTRIBUTES DLLEXPORT :: memcof
   
    ! Argument Types
    integer(C_INT), value :: n, m
    real(C_DOUBLE), intent(out) :: xms
    real(C_DOUBLE), dimension(*), intent(in) :: tseries
    real(C_DOUBLE), dimension(*), intent(out) :: d

    ! Local variables - Must match C_DOUBLE for precision
    integer(C_INT) :: i, j, k
    real(C_DOUBLE) :: p, pneum, denom
    real(C_DOUBLE), allocatable :: wk1(:), wk2(:), wkm(:)

    allocate(wk1(n), wk2(n), wkm(m))

    ! Initial power estimate
    p = 0.0_C_DOUBLE
    do j = 1, n
        p = p + tseries(j)**2
    end do
    xms = p / n

    ! Initialize forward/backward prediction arrays
    wk1(1:n) = tseries(1:n)
    wk2(1:n-1) = tseries(2:n)

    ! Burg recursion
    do k = 1, m
        pneum = 0.0_C_DOUBLE
        denom = 0.0_C_DOUBLE
       
        do j = 1, n - k
            pneum = pneum + wk1(j) * wk2(j)
            denom = denom + wk1(j)**2 + wk2(j)**2
        end do
       
        d(k) = 2.0_C_DOUBLE * pneum / denom
        xms = xms * (1.0_C_DOUBLE - d(k)**2)

        if (k > 1) then
            do i = 1, k-1
                d(i) = wkm(i) - d(k) * wkm(k-i)
            end do
        end if

        if (k == m) exit

        wkm(1:k) = d(1:k)

        do j = 1, n - k - 1
            wk1(j) = wk1(j) - wkm(k) * wk2(j)
            wk2(j) = wk2(j+1) - wkm(k) * wk1(j+1)
        end do
    end do

    deallocate(wk1, wk2, wkm)
end subroutine memcof
   
   
subroutine SymGFilter(tseries, n, lag, width, fout) BIND(C, NAME="SymGFilter")
  implicit none
  !GCC$ ATTRIBUTES DLLEXPORT :: SymGFilter
 
  integer(C_INT), value :: n, lag
  real(C_DOUBLE), intent(in) :: width
  real(C_DOUBLE), dimension(*), intent(in) :: tseries
  real(C_DOUBLE), dimension(*), intent(out) :: fout

  integer(C_INT) :: i, k, fil_len
  real(C_DOUBLE) :: sigma, norm, pi
  real(C_DOUBLE), allocatable :: kernel(:)
  pi = 3.141592653589793_C_DOUBLE

  ! Filter length = 2*lag - 1
  fil_len = 2*lag - 1

  allocate(kernel(fil_len))

  ! Compute sigma
  sigma = real(lag, kind = C_DOUBLE)/width

  ! Build Gaussian kernel centered at index lag
  do k = 1, fil_len
     kernel(k) = exp( - (real(k - lag, kind = C_DOUBLE))**2 / (2.0_C_DOUBLE * sigma**2) ) &
                 / (sqrt(2.0_C_DOUBLE * pi) * sigma)
  end do

  ! Normalize kernel
  norm = 0.0_C_DOUBLE
  do k = 1, fil_len
     norm = norm + kernel(k)
  end do
  kernel = kernel / norm

  ! Initialize output
  do i = 1, n
      fout(i) = 0.0_C_DOUBLE
  end do

  ! Convolution: valid region only
  do i = lag, n - lag + 1
     fout(i + lag - 1) = sum( kernel(:) * tseries(i - lag + 1 : i + lag - 1) )
  end do

  deallocate(kernel)
end subroutine SymGFilter
   
end module burgMemcof