chemical-equation-0.0: src/FLINT.hs

module Main where

import Common

import qualified FLINT.Utility as F
import qualified FLINT.MatrixInteger as FmpzMat
import qualified FLINT.Type as FLINT
import FLINT.Example ()

import Foreign.Ptr (Ptr)

import qualified Data.Array.Comfort.Shape as Shape
import qualified Data.Traversable as Trav
import qualified Data.Foldable as Fold
import qualified Data.Map as Map
import qualified Data.Set as Set
import Data.Set (Set)

import qualified Shell.Utility.Verbosity as Verbosity

import System.Environment (getArgs)

import Control.Monad.IO.Class (liftIO)
import Control.Monad.Trans.Cont (ContT(ContT), evalContT)



matrixFromReactantSet :: Set Reactant -> ContT r IO (Ptr FLINT.Matrix)
matrixFromReactantSet reactantSet = do
   let elements = foldMap Map.keysSet reactantSet
   mat <-
      ContT $
      F.withMatrix
         (fromIntegral $ Set.size elements)
         (fromIntegral $ Set.size reactantSet)
   liftIO $ do
      FmpzMat.zero mat
      Fold.for_ (zip [0..] $ Set.toAscList reactantSet) $ \(j,reactant) ->
         Fold.sequence_ $ flip Map.mapWithKey reactant $ \element n -> do
            let i = fromIntegral $ Shape.offset elements element
            ptr <- FmpzMat.entry mat i j
            F.integerFlintFromPrelude ptr n
   return mat


{- |
The solution consists essentially of calling 'FmpzMat.nullspace'.
All the Haskell code is needed for transferring data
between FLINT and Haskell data structures,
parsing and validating arguments and presenting the results.
-}
main :: IO ()
main = do
   let verbosity = Verbosity.normal
   (reactantSet,reactantsWithoutDuplicates)
      <- preprocessArguments verbosity =<< getArgs

   sols <- evalContT $ do
      mat <- matrixFromReactantSet reactantSet
      n <- liftIO $ FmpzMat.ncols mat
      nullBasis <- ContT $ F.withMatrix n n
      nullspaceDim <- liftIO $ FmpzMat.nullspace nullBasis mat
      let reactantIndices = Map.fromList $ zip (Set.toAscList reactantSet) [0..]
      liftIO $
         Trav.for [0 .. nullspaceDim-1] $ \j ->
            Trav.for reactantIndices $ \i ->
               F.integerPreludeFromFlint =<< FmpzMat.entry nullBasis i j

   displaySolutions verbosity reactantsWithoutDuplicates $
      map cancelIntegerVector sols