3

Intro

I'm wrapping a C numerical library with inline-c; some functions can be passed callbacks to step routines, think optimization or time integration of ODEs.

In particular in the native C use the callbacks can operate on contiguous arrays, modify them via a pointer, and return them to some opaque (distributed) datastructure.

So it's a mutable data problem, and I'd like to represent it on the Haskell side: in my understanding, in the callback we should freeze the array, work on it e.g. as a Data.Vector.Storable.Vector or with repa, thaw the result, get the foreign pointer and pass it back.

The internals: newtype Vec = Vec (Ptr Vec) deriving Storable, and the associated entry in an inline-c Context, represents the type of a pointer to an opaque C data structure and vecGetArray/vecRestoreArray produce/request the same pointer to contiguous memory and request/produce a Vec, respectively.

Q:

I noticed that, while the returned Vector is correct, when I use the resulting, modified Vec (the "side effect"), after returning from this function, it is not modified. GHC doesn't recompute it (laziness?). How do I make it recompute it? What idiomatic way is there in Haskell to work with mutable data across the FFI?


* Fixed *

see answer


Thanks!

import qualified Data.Vector.Storable as V
import qualified Data.Vector.Storable.Mutable as VM

withVecGetVectorM ::
   Vec ->                                               
   (V.Vector PetscScalar_ -> IO (V.Vector PetscScalar_)) ->   
   IO (V.Vector PetscScalar_)                           
 withVecGetVectorM v f = do 
   p <- vecGetArrayPtr v
   pf <- newForeignPtr_ p
   vImm <- V.freeze (VM.unsafeFromForeignPtr0 pf len)
   vImmOut <- f vImm        
   vMutOut <- V.thaw vImmOut
   let (fpOut, _, _) = VM.unsafeToForeignPtr vMutOut
       pOut = unsafeForeignPtrToPtr fpOut
   vecRestoreArrayPtr v pOut
   return vImmOut
where len = vecSize v

Vec.hs :

vecGetArrayPtr :: Vec -> IO (Ptr PetscScalar_)
vecGetArrayPtr v = chk1 (vecGetArrayPtr' v)

vecRestoreArrayPtr :: Vec -> Ptr PetscScalar_ -> IO ()
vecRestoreArrayPtr v ar = chk0 (vecRestoreArrayPtr' v ar)

InlineC.hs

-- PETSC_EXTERN PetscErrorCode VecGetArray(Vec,PetscScalar**);
vecGetArrayPtr' :: Vec -> IO (Ptr PetscScalar_, CInt)
vecGetArrayPtr' v = withPtr $ \p -> vga v p where
  vga v p = [C.exp|int{VecGetArray($(Vec v), $(PetscScalar** p))}|]


-- PETSC_EXTERN PetscErrorCode VecRestoreArray(Vec,PetscScalar**);
vecRestoreArrayPtr' :: Vec -> Ptr PetscScalar_ -> IO CInt
vecRestoreArrayPtr' v c = with c $ \pc -> vra v pc
  where
    vra w pc = [C.exp|int{VecRestoreArray($(Vec w), $(PetscScalar** pc))}|]

Moreover, IIUC, the code makes 2 additional copies of the vector, one at the freeze and one at thaw. , but I suspect it's inefficient. Can someone suggest improvements or simplifications?

ocramz
  • 816
  • 6
  • 18
  • The pointer inside your `Vec` type seems to be updated by `vecRestoreArrayPtr`. The error is likely in this function, so you should include it as well. – user2407038 Oct 29 '15 at 16:01
  • Thank you user2407038 , edited – ocramz Oct 29 '15 at 16:07
  • 1
    Also, you can get rid of the copies by replacing `thaw` with `unsafeThaw`, which doesn't make a copy. This should be safe, I think, since you overwrite your `Vec` with the result of the computation anyways. However, if an exception occurs, you may be left in an inconsistent state - might be better to copy it over anyways. – user2407038 Oct 29 '15 at 16:10

1 Answers1

0

Fixed:

I was making a trivial mistake re. which pointers to handle; here's the fix. This function has the same signature as the buggy one, which suggests me there has to be a typeful way to signify which pointer we are modifying. If anyone has suggestions regarding this, do not hesitate to get in touch.

vecRestoreVector :: Vec -> V.Vector PetscScalar_ -> IO ()
vecRestoreVector v w = do
  p <- vecGetArrayPtr v
  pf <- newForeignPtr_ p
  V.copy (VM.unsafeFromForeignPtr0 pf len) w
  vecRestoreArrayPtr v p
    where
     len = vecSize v
Community
  • 1
  • 1
ocramz
  • 816
  • 6
  • 18