{-# LINE 1 "src/Data/Number/Flint/Support/Mpf/Mat/FFI.hsc" #-}
{-|
module      :  Data.Number.Flint.Support.Mpf.Mat.FFI
copyright   :  (c) 2022 Hartmut Monien
license     :  GNU GPL, version 2 or above (see LICENSE)
maintainer  :  hmonien@uni-bonn.de
-}
module Data.Number.Flint.Support.Mpf.Mat.FFI (
  -- * Matrices of MPF floating-point numbers
    MpfMat (..)
  , CMpfMat (..)
  , newMpfMat
  , withMpfMat
  , withNewMpfMat
  -- * Memory management
  , mpf_mat_init
  , mpf_mat_clear
  -- * Basic assignment and manipulation
  , mpf_mat_set
  , mpf_mat_swap
  , mpf_mat_swap_entrywise
  , mpf_mat_entry
  , mpf_mat_zero
  , mpf_mat_one
  -- * Random matrix generation
  , mpf_mat_randtest
  -- * Input and output
  , mpf_mat_print
  -- * Comparison
  , mpf_mat_equal
  , mpf_mat_approx_equal
  , mpf_mat_is_zero
  , mpf_mat_is_empty
  , mpf_mat_is_square
  -- * Matrix multiplication
  , mpf_mat_mul
  -- * Gram-Schmidt Orthogonalisation and QR Decomposition
  , mpf_mat_gso
  , mpf_mat_qr
) where 

-- matrices of MPF floating-point numbers --------------------------------------

import Foreign.C.Types
import Foreign.Ptr
import Foreign.ForeignPtr
import Foreign.Storable
import Foreign.Marshal.Array

import Data.Number.Flint.Flint




-- mpf_mat_t -------------------------------------------------------------------

data MpfMat = MpfMat {-# UNPACK #-} !(ForeignPtr CMpfMat)
data CMpfMat = CMpfMat (Ptr CMpf) CLong CLong (Ptr (Ptr CMpf))

instance Storable CMpfMat where
  {-# INLINE sizeOf #-}
  sizeOf :: CMpfMat -> Int
sizeOf CMpfMat
_ = (Int
40)
{-# LINE 62 "src/Data/Number/Flint/Support/Mpf/Mat/FFI.hsc" #-}
  {-# INLINE alignment #-}
  alignment :: CMpfMat -> Int
alignment CMpfMat
_ = Int
8
{-# LINE 64 "src/Data/Number/Flint/Support/Mpf/Mat/FFI.hsc" #-}
  peek ptr = CMpfMat
    <$> (\hsc_ptr -> peekByteOff hsc_ptr 0) ptr
{-# LINE 66 "src/Data/Number/Flint/Support/Mpf/Mat/FFI.hsc" #-}
    <*> (\hsc_ptr -> peekByteOff hsc_ptr 8) ptr
{-# LINE 67 "src/Data/Number/Flint/Support/Mpf/Mat/FFI.hsc" #-}
    <*> (\hsc_ptr -> peekByteOff hsc_ptr 16) ptr
{-# LINE 68 "src/Data/Number/Flint/Support/Mpf/Mat/FFI.hsc" #-}
    <*> (\hsc_ptr -> peekByteOff hsc_ptr 32) ptr
{-# LINE 69 "src/Data/Number/Flint/Support/Mpf/Mat/FFI.hsc" #-}
  poke = error "CMpfMat.poke: Not defined."

newMpfMat :: CLong -> CLong -> CFBitCnt -> IO MpfMat
newMpfMat CLong
rows CLong
cols CFBitCnt
prec = do
  ForeignPtr CMpfMat
x <- IO (ForeignPtr CMpfMat)
forall a. Storable a => IO (ForeignPtr a)
mallocForeignPtr
  ForeignPtr CMpfMat -> (Ptr CMpfMat -> IO ()) -> IO ()
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr CMpfMat
x ((Ptr CMpfMat -> IO ()) -> IO ())
-> (Ptr CMpfMat -> IO ()) -> IO ()
forall a b. (a -> b) -> a -> b
$ \Ptr CMpfMat
x -> Ptr CMpfMat -> CLong -> CLong -> CFBitCnt -> IO ()
mpf_mat_init Ptr CMpfMat
x CLong
rows CLong
cols CFBitCnt
prec
  FinalizerPtr CMpfMat -> ForeignPtr CMpfMat -> IO ()
forall a. FinalizerPtr a -> ForeignPtr a -> IO ()
addForeignPtrFinalizer FinalizerPtr CMpfMat
p_mpf_mat_clear ForeignPtr CMpfMat
x
  MpfMat -> IO MpfMat
forall a. a -> IO a
forall (m :: * -> *) a. Monad m => a -> m a
return (MpfMat -> IO MpfMat) -> MpfMat -> IO MpfMat
forall a b. (a -> b) -> a -> b
$ ForeignPtr CMpfMat -> MpfMat
MpfMat ForeignPtr CMpfMat
x

{-# INLINE withMpfMat #-}
withMpfMat :: MpfMat -> (Ptr CMpfMat -> IO a) -> IO (MpfMat, a)
withMpfMat (MpfMat ForeignPtr CMpfMat
x) Ptr CMpfMat -> IO a
f = do
  ForeignPtr CMpfMat
-> (Ptr CMpfMat -> IO (MpfMat, a)) -> IO (MpfMat, a)
forall a b. ForeignPtr a -> (Ptr a -> IO b) -> IO b
withForeignPtr ForeignPtr CMpfMat
x ((Ptr CMpfMat -> IO (MpfMat, a)) -> IO (MpfMat, a))
-> (Ptr CMpfMat -> IO (MpfMat, a)) -> IO (MpfMat, a)
forall a b. (a -> b) -> a -> b
$ \Ptr CMpfMat
px -> Ptr CMpfMat -> IO a
f Ptr CMpfMat
px IO a -> (a -> IO (MpfMat, a)) -> IO (MpfMat, a)
forall a b. IO a -> (a -> IO b) -> IO b
forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
>>= (MpfMat, a) -> IO (MpfMat, a)
forall a. a -> IO a
forall (m :: * -> *) a. Monad m => a -> m a
return ((MpfMat, a) -> IO (MpfMat, a))
-> (a -> (MpfMat, a)) -> a -> IO (MpfMat, a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (ForeignPtr CMpfMat -> MpfMat
MpfMat ForeignPtr CMpfMat
x,)

withNewMpfMat :: CLong
-> CLong -> CFBitCnt -> (Ptr CMpfMat -> IO a) -> IO (MpfMat, a)
withNewMpfMat CLong
rows CLong
cols CFBitCnt
prec Ptr CMpfMat -> IO a
f = do
  MpfMat
x <- CLong -> CLong -> CFBitCnt -> IO MpfMat
newMpfMat CLong
rows CLong
cols CFBitCnt
prec
  MpfMat -> (Ptr CMpfMat -> IO a) -> IO (MpfMat, a)
forall {a}. MpfMat -> (Ptr CMpfMat -> IO a) -> IO (MpfMat, a)
withMpfMat MpfMat
x Ptr CMpfMat -> IO a
f
  
-- Memory management -----------------------------------------------------------

-- | /mpf_mat_init/ /mat/ /rows/ /cols/ /prec/ 
-- 
-- Initialises a matrix with the given number of rows and columns and the
-- given precision for use. The precision is at least the precision of the
-- entries.
foreign import ccall "mpf_mat.h mpf_mat_init"
  mpf_mat_init :: Ptr CMpfMat -> CLong -> CLong -> CFBitCnt -> IO ()

-- | /mpf_mat_clear/ /mat/ 
-- 
-- Clears the given matrix.
foreign import ccall "mpf_mat.h mpf_mat_clear"
  mpf_mat_clear :: Ptr CMpfMat -> IO ()

foreign import ccall "mpf_mat.h &mpf_mat_clear"
  p_mpf_mat_clear :: FunPtr (Ptr CMpfMat -> IO ())

-- Basic assignment and manipulation -------------------------------------------

-- | /mpf_mat_set/ /mat1/ /mat2/ 
-- 
-- Sets @mat1@ to a copy of @mat2@. The dimensions of @mat1@ and @mat2@
-- must be the same.
foreign import ccall "mpf_mat.h mpf_mat_set"
  mpf_mat_set :: Ptr CMpfMat -> Ptr CMpfMat -> IO ()

-- | /mpf_mat_swap/ /mat1/ /mat2/ 
-- 
-- Swaps two matrices. The dimensions of @mat1@ and @mat2@ are allowed to
-- be different.
foreign import ccall "mpf_mat.h mpf_mat_swap"
  mpf_mat_swap :: Ptr CMpfMat -> Ptr CMpfMat -> IO ()

-- | /mpf_mat_swap_entrywise/ /mat1/ /mat2/ 
-- 
-- Swaps two matrices by swapping the individual entries rather than
-- swapping the contents of the structs.
foreign import ccall "mpf_mat.h mpf_mat_swap_entrywise"
  mpf_mat_swap_entrywise :: Ptr CMpfMat -> Ptr CMpfMat -> IO ()

-- | /mpf_mat_entry/ /mat/ /i/ /j/ 
-- 
-- Returns a reference to the entry of @mat@ at row \(i\) and column \(j\).
-- Both \(i\) and \(j\) must not exceed the dimensions of the matrix. The
-- return value can be used to either retrieve or set the given entry.
mpf_mat_entry :: Ptr CMpfMat -> CLong -> CLong -> IO (Ptr CMpf)
mpf_mat_entry :: Ptr CMpfMat -> CLong -> CLong -> IO (Ptr CMpf)
mpf_mat_entry Ptr CMpfMat
mat CLong
i CLong
j = do
  CMpfMat Ptr CMpf
entries CLong
r CLong
c Ptr (Ptr CMpf)
rows <- Ptr CMpfMat -> IO CMpfMat
forall a. Storable a => Ptr a -> IO a
peek Ptr CMpfMat
mat
  Ptr CMpf -> IO (Ptr CMpf)
forall a. a -> IO a
forall (m :: * -> *) a. Monad m => a -> m a
return (Ptr CMpf -> IO (Ptr CMpf)) -> Ptr CMpf -> IO (Ptr CMpf)
forall a b. (a -> b) -> a -> b
$ Ptr CMpf
entries Ptr CMpf -> Int -> Ptr CMpf
forall a. Storable a => Ptr a -> Int -> Ptr a
`advancePtr` (CLong -> Int
forall a b. (Integral a, Num b) => a -> b
fromIntegral (CLong
iCLong -> CLong -> CLong
forall a. Num a => a -> a -> a
*CLong
c CLong -> CLong -> CLong
forall a. Num a => a -> a -> a
+ CLong
j))

-- | /mpf_mat_zero/ /mat/ 
-- 
-- Sets all entries of @mat@ to 0.
foreign import ccall "mpf_mat.h mpf_mat_zero"
  mpf_mat_zero :: Ptr CMpfMat -> IO ()

-- | /mpf_mat_one/ /mat/ 
-- 
-- Sets @mat@ to the unit matrix, having ones on the main diagonal and
-- zeroes elsewhere. If @mat@ is nonsquare, it is set to the truncation of
-- a unit matrix.
foreign import ccall "mpf_mat.h mpf_mat_one"
  mpf_mat_one :: Ptr CMpfMat -> IO ()

-- Random matrix generation ----------------------------------------------------

-- | /mpf_mat_randtest/ /mat/ /state/ /bits/ 
-- 
-- Sets the entries of @mat@ to random numbers in the interval \([0, 1)\)
-- with @bits@ significant bits in the mantissa or less if their precision
-- is smaller.
foreign import ccall "mpf_mat.h mpf_mat_randtest"
  mpf_mat_randtest :: Ptr CMpfMat -> Ptr CFRandState -> CFBitCnt -> IO ()

-- Input and output ------------------------------------------------------------

-- | /mpf_mat_print/ /mat/ 
-- 
-- Prints the given matrix to the stream @stdout@.
foreign import ccall "mpf_mat.h mpf_mat_print"
  mpf_mat_print :: Ptr CMpfMat -> IO ()

-- Comparison ------------------------------------------------------------------

-- | /mpf_mat_equal/ /mat1/ /mat2/ 
-- 
-- Returns a non-zero value if @mat1@ and @mat2@ have the same dimensions
-- and entries, and zero otherwise.
foreign import ccall "mpf_mat.h mpf_mat_equal"
  mpf_mat_equal :: Ptr CMpfMat -> Ptr CMpfMat -> IO CInt

-- | /mpf_mat_approx_equal/ /mat1/ /mat2/ /bits/ 
-- 
-- Returns a non-zero value if @mat1@ and @mat2@ have the same dimensions
-- and the first @bits@ bits of their entries are equal, and zero
-- otherwise.
foreign import ccall "mpf_mat.h mpf_mat_approx_equal"
  mpf_mat_approx_equal :: Ptr CMpfMat -> Ptr CMpfMat -> CFBitCnt -> IO CInt

-- | /mpf_mat_is_zero/ /mat/ 
-- 
-- Returns a non-zero value if all entries @mat@ are zero, and otherwise
-- returns zero.
foreign import ccall "mpf_mat.h mpf_mat_is_zero"
  mpf_mat_is_zero :: Ptr CMpfMat -> IO CInt

-- | /mpf_mat_is_empty/ /mat/ 
-- 
-- Returns a non-zero value if the number of rows or the number of columns
-- in @mat@ is zero, and otherwise returns zero.
foreign import ccall "mpf_mat.h mpf_mat_is_empty"
  mpf_mat_is_empty :: Ptr CMpfMat -> IO CInt

-- | /mpf_mat_is_square/ /mat/ 
-- 
-- Returns a non-zero value if the number of rows is equal to the number of
-- columns in @mat@, and otherwise returns zero.
foreign import ccall "mpf_mat.h mpf_mat_is_square"
  mpf_mat_is_square :: Ptr CMpfMat -> IO CInt

-- Matrix multiplication -------------------------------------------------------

-- | /mpf_mat_mul/ /C/ /A/ /B/ 
-- 
-- Sets @C@ to the matrix product \(C = A B\). The matrices must have
-- compatible dimensions for matrix multiplication (an exception is raised
-- otherwise). Aliasing is allowed.
foreign import ccall "mpf_mat.h mpf_mat_mul"
  mpf_mat_mul :: Ptr CMpfMat -> Ptr CMpfMat -> Ptr CMpfMat -> IO ()

-- Gram-Schmidt Orthogonalisation and QR Decomposition -------------------------

-- | /mpf_mat_gso/ /B/ /A/ 
-- 
-- Takes a subset of \(R^m\) \(S = {a_1, a_2, \ldots ,a_n}\) (as the
-- columns of a \(m x n\) matrix @A@) and generates an orthonormal set
-- \(S' = {b_1, b_2, \ldots ,b_n}\) (as the columns of the \(m x n\) matrix
-- @B@) that spans the same subspace of \(R^m\) as \(S\).
-- 
-- This uses an algorithm of Schwarz-Rutishauser. See pp. 9 of
-- <https://people.inf.ethz.ch/gander/papers/qrneu.pdf>
foreign import ccall "mpf_mat.h mpf_mat_gso"
  mpf_mat_gso :: Ptr CMpfMat -> Ptr CMpfMat -> IO ()

-- | /mpf_mat_qr/ /Q/ /R/ /A/ 
-- 
-- Computes the \(QR\) decomposition of a matrix @A@ using the Gram-Schmidt
-- process. (Sets @Q@ and @R@ such that \(A = QR\) where @R@ is an upper
-- triangular matrix and @Q@ is an orthogonal matrix.)
-- 
-- This uses an algorithm of Schwarz-Rutishauser. See pp. 9 of
-- <https://people.inf.ethz.ch/gander/papers/qrneu.pdf>
foreign import ccall "mpf_mat.h mpf_mat_qr"
  mpf_mat_qr :: Ptr CMpfMat -> Ptr CMpfMat -> Ptr CMpfMat -> IO ()