WCV example

copilot-verifier/copilot-verifier.cabal

@@ -73,6 +73,7 @@ library copilot-verifier-examples
     Copilot.Verifier.Examples.Heater
     Copilot.Verifier.Examples.Structs
     Copilot.Verifier.Examples.Voting
+    Copilot.Verifier.Examples.WCV
 
 executable verify-examples
   import: bldflags

copilot-verifier/examples/Copilot/Verifier/Examples.hs

@@ -15,6 +15,7 @@ import qualified Copilot.Verifier.Examples.Engine  as Engine
 import qualified Copilot.Verifier.Examples.Heater  as Heater
 import qualified Copilot.Verifier.Examples.Structs as Structs
 import qualified Copilot.Verifier.Examples.Voting  as Voting
+import qualified Copilot.Verifier.Examples.WCV     as WCV
 
 allExamples :: Map (CI Text) (IO ())
 allExamples = Map.fromList
@@ -26,6 +27,7 @@ allExamples = Map.fromList
     , example "Heater" Heater.main
     , example "Structs" Structs.main
     , example "Voting" Voting.main
+    , example "WCV"    WCV.main
     ]
   where
     example :: Text -> IO () -> (CI Text, IO ())

copilot-verifier/examples/Copilot/Verifier/Examples/WCV.hs

@@ -0,0 +1,192 @@
+-- | This example shows an implementation of the Well-Clear Violation
+-- algorithm, it follows the implementation described in 'Analysis of
+-- Well-Clear Bounday Models for the Integration of UAS in the NAS',
+-- https://ntrs.nasa.gov/citations/20140010078.
+
+{-# LANGUAGE DataKinds #-}
+{-# LANGUAGE RebindableSyntax #-}
+
+module Copilot.Verifier.Examples.WCV where
+
+import Language.Copilot
+import qualified Copilot.Theorem.What4 as CT
+import Copilot.Compile.C99
+import Copilot.Verifier (verify)
+
+import qualified Prelude as P
+import Data.Foldable (forM_)
+import qualified Control.Monad as Monad
+
+
+-- | `dthr` is the horizontal distance threshold.
+dthr :: Stream Double
+dthr = extern "dthr" Nothing
+
+-- | `tthr` is the horizontal time threshold.
+tthr :: Stream Double
+tthr = extern "tthr" Nothing
+
+-- | `zthr` is the vertical distance / altitude threshold.
+zthr :: Stream Double
+zthr = extern "zthr" Nothing
+
+-- | `tcoathr` is the vertical time threshold.
+tcoathr :: Stream Double
+tcoathr = extern "tcoathr" Nothing
+
+type Vect2 = (Stream Double, Stream Double)
+
+
+--------------------------------
+-- External streams for relative position and velocity.
+--------------------------------
+
+-- | The relative x velocity between ownship and the intruder.
+vx :: Stream Double
+vx = extern "relative_velocity_x" Nothing
+
+-- | The relative y velocity between ownship and the intruder.
+vy :: Stream Double
+vy = extern "relative_velocity_y" Nothing
+
+-- | The relative z velocity between ownship and the intruder.
+vz :: Stream Double
+vz = extern "relative_velocity_z" Nothing
+
+-- | The relative velocity as a 2D vector.
+v :: (Stream Double, Stream Double)
+v = (vx, vy)
+
+
+-- | The relative x position between ownship and the intruder.
+sx :: Stream Double
+sx = extern "relative_position_x" Nothing
+
+-- | The relative y position between ownship and the intruder.
+sy :: Stream Double
+sy = extern "relative_position_y" Nothing
+
+-- | The relative z position between ownship and the intruder.
+sz :: Stream Double
+sz = extern "relative_position_z" Nothing
+
+-- | The relative position as a 2D vector.
+s :: (Stream Double, Stream Double)
+s = (sx, sy)
+
+
+------------------
+-- The following section contains basic libraries for working with vectors.
+------------------
+
+-- | Multiply two Vectors.
+(|*|) :: Vect2 -> Vect2 -> Stream Double
+(|*|) (x1, y1) (x2, y2) = (x1 * x2) + (y1 * y2)
+
+-- | Calculate the square of a vector.
+sq :: Vect2 -> Stream Double
+sq x = x |*| x
+
+-- | Calculate the length of a vector.
+norm :: Vect2 -> Stream Double
+norm = sqrt . sq
+
+-- | Calculate the determinant of two vectors.
+det :: Vect2 -> Vect2 -> Stream Double
+det (x1, y1) (x2, y2) = (x1 * y2) - (x2 * y1)
+
+-- | Compare two vectors, taking into account the small error that is
+-- introduced by the usage of `Double`s.
+(~=) :: Stream Double -> Stream Double -> Stream Bool
+a ~= b = (abs (a - b)) < 0.001
+
+-- | Negate a vector.
+neg :: Vect2 -> Vect2
+neg (x, y) = (negate x, negate y)
+
+
+--------------------
+-- From here on the algorithm, as described by the paper mentioned on the top
+-- of this file, is implemented. Please refer to the paper for details.
+--------------------
+
+tau :: Vect2 -> Vect2 -> Stream Double
+tau s v = if s |*| v < 0
+            then (-(sq s)) / (s |*| v)
+            else -1
+
+tcpa :: Vect2 -> Vect2 -> Stream Double
+tcpa s v@(vx, vy) = if vx ~= 0 && vy ~= 0
+                      then 0
+                      else -(s |*| v)/(sq v)
+
+taumod :: Vect2 -> Vect2 -> Stream Double
+taumod s v = if s |*| v < 0
+               then (dthr * dthr - (sq s))/(s |*| v)
+               else -1
+
+tep :: Vect2 -> Vect2 -> Stream Double
+tep s v = if (s |*| v < 0) && ((delta s v dthr) >= 0)
+            then theta s v dthr (-1)
+            else -1
+
+delta :: Vect2 -> Vect2 -> Stream Double -> Stream Double
+delta s v d = (d*d) * (sq v) - ((det s v)*(det s v))
+-- Here the formula says : (s . orth v)^2 which is the same as det(s,v)^2
+
+theta :: Vect2 -> Vect2 -> Stream Double -> Stream Double -> Stream Double
+theta s v d e = (-(s |*| v) + e * (sqrt $ delta s v d)) / (sq v)
+
+
+tcoa :: Stream Double -> Stream Double -> Stream Double
+tcoa sz vz = if (sz * vz) < 0
+               then (-sz) / vz
+               else -1
+
+dcpa :: Vect2 -> Vect2 -> Stream Double
+dcpa s@(sx, sy) v@(vx, vy) = norm (sx + (tcpa s v) * vx, sy + (tcpa s v) * vy)
+
+
+--------------------------
+-- Well clear Violation --
+--------------------------
+
+-- | Determines if the well clear property is violated or not.
+wcv :: (Vect2 -> Vect2 -> Stream Double) ->
+       Vect2 -> Stream Double ->
+       Vect2 -> Stream Double ->
+       Stream Bool
+wcv tvar s sz v vz = (horizontalWCV tvar s v) && (verticalWCV sz vz)
+
+verticalWCV :: Stream Double -> Stream Double -> Stream Bool
+verticalWCV sz vz =
+  ((abs $ sz) <= zthr) ||
+  (0 <= (tcoa sz vz) && (tcoa sz vz) <= tcoathr)
+
+horizontalWCV :: (Vect2 -> Vect2 -> Stream Double) -> Vect2 -> Vect2 -> Stream Bool
+horizontalWCV tvar s v =
+  (norm s <= dthr) ||
+  (((dcpa s v) <= dthr) && (0 <= (tvar s v)) && ((tvar s v) <= tthr))
+
+spec = do
+  Monad.void $ prop "1a" (forall $ (tau s v) ~= (tau (neg s) (neg v)))
+  -- Monad.void $ prop "3d" (forall $ (wcv tep s sz v vz)    == (wcv tep (neg s) (-sz) (neg v) (-vz)))
+  trigger "well_clear_violation" (wcv tep s sz v vz) []
+
+main :: IO ()
+main = do
+  spec' <- reify spec
+  verify mkDefaultCSettings [] "wcv" spec'
+
+{-
+  -- Use Z3 to prove the properties.
+  results <- CT.prove CT.Z3 spec'
+
+  -- Print the results.
+  forM_ results $ \(nm, res) -> do
+    putStr $ nm <> ": "
+    case res of
+      CT.Valid -> putStrLn "valid"
+      CT.Invalid -> putStrLn "invalid"
+      CT.Unknown -> putStrLn "unknown"
+-}

copilot-verifier/src/Copilot/Verifier.hs

@@ -47,7 +47,8 @@ import Lang.Crucible.Backend
   , pushAssumptionFrame, popUntilAssumptionFrame
   , getProofObligations, clearProofObligations
   , LabeledPred(..), addDurableProofObligation
-  , addAssumption, CrucibleAssumption(..)
+  , addAssumption, CrucibleAssumption(..), ppAbortExecReason
+  -- , ProofObligations, proofGoal, goalsToList, labeledPredMsg
   )
 import Lang.Crucible.Backend.Simple (newSimpleBackend)
 import Lang.Crucible.CFG.Core (AnyCFG(..), cfgArgTypes, cfgReturnType)
@@ -58,10 +59,11 @@ import Lang.Crucible.Simulator
   , defaultAbortHandler, runOverrideSim, partialValue, gpValue
   , GlobalVar, executeCrucible, OverrideSim, regValue
   , readGlobal, writeGlobal, callCFG, emptyRegMap, RegEntry(..)
+  , AbortedResult(..)
   )
 import Lang.Crucible.Simulator.GlobalState ( insertGlobal )
 import Lang.Crucible.Simulator.RegValue (RegValue, RegValue'(..))
-import Lang.Crucible.Simulator.SimError (SimError(..), SimErrorReason(..))
+import Lang.Crucible.Simulator.SimError (SimError(..), SimErrorReason(..)) -- ppSimError
 import Lang.Crucible.Types
   ( TypeRepr(..), (:~:)(..), KnownRepr(..), BoolType )
 
@@ -96,7 +98,7 @@ import Lang.Crucible.LLVM.TypeContext (TypeContext, llvmDataLayout)
 import Crux (defaultOutputConfig)
 import Crux.Config (cfgJoin, Config(..))
 import Crux.Config.Load (fromFile, fromEnv)
-import Crux.Config.Common (cruxOptions, CruxOptions(..), postprocessOptions, )
+import Crux.Config.Common (cruxOptions, CruxOptions(..), postprocessOptions, outputOptions )
 import Crux.Goal (proveGoalsOffline, provedGoalsTree)
 import Crux.Log
   ( cruxLogMessageToSayWhat, withCruxLogMessage, outputHandle
@@ -121,6 +123,10 @@ import What4.Interface
   , truePred, falsePred, eqPred, andPred, backendPred
   )
 import What4.Expr.Builder (FloatModeRepr(..), ExprBuilder, BoolExpr, startCaching)
+import What4.InterpretedFloatingPoint
+  ( FloatInfoRepr(..), IsInterpretedFloatExprBuilder(..)
+  , SingleFloat, DoubleFloat
+  )
 import What4.ProgramLoc (ProgramLoc, mkProgramLoc, Position(..))
 import What4.Solver.Adapter (SolverAdapter(..))
 import What4.Solver.Z3 (z3Adapter)
@@ -175,6 +181,8 @@ verifyBitcode csettings properties spec cruxOpts llvmOpts bcFile =
      startCaching sym
      bbMapRef <- newIORef mempty
      let ?recordLLVMAnnotation = \an bb -> modifyIORef bbMapRef (Map.insert an bb)
+     ocfg <- defaultOutputConfig cruxLLVMLoggingToSayWhat
+     let ?outputConfig = ocfg (Just (outputOptions cruxOpts))
 
      let adapters = [z3Adapter] -- TODO? configurable
      extendConfig (solver_adapter_config_options z3Adapter) (getConfiguration sym)
@@ -377,7 +385,12 @@ executeStep csettings simctx memVar mem llvmmod modTrans triggerGlobals triggerO
      res <- executeCrucible [] initSt
      case res of
        FinishedResult _ pr -> return (pr^.partialValue.gpValue.to regValue)
-       AbortedResult{} -> fail "simulation aborted!"
+       AbortedResult _ (AbortedExec rsn _) ->
+         fail $ unlines $ ["Simulation aborted!", show (ppAbortExecReason rsn)]
+       AbortedResult _ (AbortedExit ec) ->
+         fail $ unlines $ ["Simulation called exit!", show ec]
+       AbortedResult _ (AbortedBranch{}) ->
+         fail $ unlines $ ["Simulation aborted for multiple reasons."]
        TimeoutResult{} -> fail "simulation timed out!"
  where
   setupTrigger gs (_,gv,_) = insertGlobal gv (falsePred sym) gs
@@ -591,10 +604,8 @@ copilotExprToRegValue sym = loop
     loop (CW4.XWord64 x) (LLVMPointerRepr w) | Just Refl <- testEquality w (knownNat @64) =
       llvmPointer_bv sym x
 
---    loop (CW4.XFloat x)  (FloatRepr SingleFloatRepr) v =
---      iFloatEq sym x v
---    loop (CW4.XDouble x) (FloatRepr DoubleFloatRepr) v =
---      iFloatEq sym x v
+    loop (CW4.XFloat x)  (FloatRepr SingleFloatRepr) = return x
+    loop (CW4.XDouble x) (FloatRepr DoubleFloatRepr) = return x
 
     loop CW4.XEmptyArray (VectorRepr _tpr) =
       pure V.empty
@@ -646,11 +657,8 @@ computeEqualVals sym mem = loop
     loop Word64 (CW4.XWord64 x) (LLVMPointerRepr w) v | Just Refl <- testEquality w (knownNat @64) =
       bvEq sym x =<< projectLLVM_bv sym v
 
-
---    loop Float (CW4.XFloat x)  (FloatRepr SingleFloatRepr) v =
---      iFloatEq sym x v
---    loop Double (CW4.XDouble x) (FloatRepr DoubleFloatRepr) v =
---      iFloatEq sym x v
+    loop Float (CW4.XFloat x)  (FloatRepr SingleFloatRepr) v = iFloatEq @_ @SingleFloat sym x v
+    loop Double (CW4.XDouble x) (FloatRepr DoubleFloatRepr) v = iFloatEq @_ @DoubleFloat sym x v
 
     loop (Array _ctp) CW4.XEmptyArray (VectorRepr _tpr) vs =
       pure $ backendPred sym $ V.null vs
@@ -866,9 +874,17 @@ proveObls cruxOpts adapters bbMapRef simctx =
      obls <- getProofObligations sym
      clearProofObligations sym
 
+--     mapM_ (print . ppSimError) (summarizeObls sym obls)
+
      results <- proveGoalsOffline adapters cruxOpts simctx (explainFailure sym bbMapRef) obls
      presentResults sym results
 
+{-
+summarizeObls :: sym -> ProofObligations sym -> [SimError]
+summarizeObls _ Nothing = []
+summarizeObls _ (Just obls) = map (view labeledPredMsg . proofGoal) (goalsToList obls)
+-}
+
 presentResults ::
   Logs msgs =>
   IsSymInterface sym =>