Haskellは素晴らしい

何故か、Vertigoのソースを眺めていたら(読めない)、IOExtsをimportしていたので、調べて見たら、unsafePerformIO :: IO a -> aとかいう関数が。一旦、モナド界に行ったら、現世には帰ってこれないはずなんじゃ･･･。Maybeモナドにも、fromJustとかあるし。むー

import System.IO.Unsafe
main=
 let c = (unsafePerformIO getChar)
  in if (c=='a') then (putChar 'A') else (putChar c)

多項式計算。多項式の正規化(normalize)が全て。もっと簡単に書けると思って始めたら、すごい長くなった。もっと賢いやり方があるはず

module Polynomial where

data Polynomial a = 
 Const a | Var String |PLUS (Polynomial a) (Polynomial a)|MULT (Polynomial a) (Polynomial a)


instance Show a => Show (Polynomial a) where
  show (Const n) = show n
  show (Var str) = show str
  show (PLUS x y) = "(+ " ++ " " ++ (show x) ++ " " ++ (show y) ++ ")"
  show (MULT x y) = "(* " ++ " " ++ (show x) ++ " " ++ (show y) ++ ")"

instance Num a => Eq (Polynomial a) where
   p == q = (term_eq (normalize p) (normalize q))

instance Num a => Num (Polynomial a) where
   p+q = (normalize (PLUS p q))
   p*q = (normalize (MULT p q))
   p-q = (normalize (PLUS p (MULT (Const (negate 1)) q)))
   fromInteger n = (Const (fromInteger n))

--項としての比較。x+4と4+xは多項式としては等しいが、項としては異なる
--多項式としての比較は、Eqのinstanceとして
term_eq :: (Eq a) => (Polynomial a) -> (Polynomial a) -> Bool
term_eq p q = case (p,q) of
  ((Const x) , (Const y))   ->  (x == y)
  ((Var x)  , (Var y)  )   -> (x == y)
  ((PLUS x1 x2) , (PLUS y1 y2))  ->  (term_eq x1 y1) && (term_eq x2 y2)
  ((MULT x1 x2) , (MULT y1 y2))  ->  (term_eq x1 y1) && (term_eq x2 y2)
  ( _  , _) -> False


--単項式の和に展開
expand :: (Num a) => (Polynomial a) -> (Polynomial a)
expand p = 
 let next = (expand1 p)
  in if (term_eq p next) then next else (expand next)
 where
  expand1 p = case p of
    (Const n) -> (Const n)
    (Var v) -> (Var v)
    (PLUS x y) ->  (PLUS (expand1 x) (expand1 y))
    (MULT x (PLUS y z)) -> (expand1 (PLUS (MULT x y) (MULT x z)))
    (MULT (PLUS x y) z) -> (expand1 (PLUS (MULT x z) (MULT y z)))
    (MULT x y) -> (MULT (expand1 x) (expand1 y))
  
--単項式を正規化(引数は単項式と仮定されている。頭悪い)
----単項式をリスト化(e.g. x*3*y->["x" , 3 ,"y"]
flatten_monomial :: (Num a) => (Polynomial a) -> [(Polynomial a)]
flatten_monomial p = case p of
  (Const n) -> [(Const n)]
  (Var x) -> [(Var x)]
  (MULT x y) -> (flatten_monomial x) ++ (flatten_monomial y)

----単項式p,qの比較(p<q)。正規化してあること前提。sort_monomialでも使う
gt_monomial :: (Polynomial a) -> (Polynomial a) -> Bool
gt_monomial p q = case (p,q) of
  ((Const n) , (Const m)) -> False
  ((Const n) , p) -> False
  (p , (Const n)) -> True
  ((Var x) , (Var y)) -> (x < y)
  ((Var x) , (MULT (Var z) w)) -> (x > z)
  ((MULT (Var z) w) , (Var x)) -> (not (x < z))
  ((MULT (Var x) y) , (MULT (Var z) w)) -> (if (x == z) then (gt_monomial y w) else (x < z))
  ((MULT (Const n) x) , y) -> (gt_monomial x y)
  (x , (MULT (Const n) y)) -> (gt_monomial x y)

----リスト化した単項式をソート
sort_monomial :: (Num a) => [(Polynomial a)] -> [(Polynomial a)]
sort_monomial [] = []
sort_monomial (x:xs) 
  = (sort_monomial [y|y<-xs , (gt_monomial x y)]) 
     ++ [x] ++ (sort_monomial [y|y<-xs , (not (gt_monomial x y))])

----平坦化した単項式を戻す
list_to_monomial :: (Num a) => [(Polynomial a)] -> (Polynomial a)
list_to_monomial p = case p of
  [x] -> x
  ((Const 0):xs) -> (Const 0)
  ((Const 1):xs) -> (list_to_monomial (reverse xs))
  ((Const n):xs) -> (MULT (Const n) (list_to_monomial (reverse xs)))
  (x:xs) -> (MULT x (list_to_monomial xs))

----定数部分を計算してしまう
reduce_monomial :: (Num a) => [(Polynomial a)] -> [(Polynomial a)]
reduce_monomial v = case v of
   ((Const n):((Const m):xs)) -> (reduce_monomial ((Const (n*m)) : xs))
   ((Var x):xs) -> ((Const 1):((Var x):xs))
   x ->x

----Normalizing Monomial
normalize_monomial :: (Num a) => (Polynomial a) -> (Polynomial a)
normalize_monomial p 
       = list_to_monomial (reduce_monomial (sort_monomial (flatten_monomial p)))


--多項式の正規化(頭悪い)
----多項式を(正規化した)単項式のリストに分解。expandしてあること前提
flatten_polynomial :: (Num a) => (Polynomial a) -> [(Polynomial a)]
flatten_polynomial p = case p of
  (PLUS x y) -> (flatten_polynomial x) ++ (flatten_polynomial y)
  x -> [(normalize_monomial x)]

----単項式の足し算
add_monomial :: (Num a) => (Polynomial a) -> (Polynomial a) -> (Polynomial a)
add_monomial p q = case (p,q) of
  ((Const n) , (Const m)) -> (Const (n+m))
  ((MULT (Const n) x) , (MULT (Const m) y)) 
    -> (if (eq_monomial x y)
	then (MULT (Const (n+m)) x)
	else (PLUS (MULT (Const n) x)  (MULT (Const m) y)))
  (x , (MULT (Const n) y))
    -> (if (eq_monomial x y) 
        then (MULT (Const (n+1)) x)
        else (PLUS x (MULT (Const n) y)))
  ((MULT (Const n) x) , y)
    -> (if (eq_monomial x y) 
        then (MULT (Const (n+1)) x)
        else (PLUS (MULT (Const n) x) y))
  (x , y) -> (PLUS x y)
 where
   eq_monomial x y = (not (gt_monomial x y)) && (not (gt_monomial y x))

----単項式リストを多項式に戻す
list_to_polynomial :: (Num a) => [(Polynomial a)] -> (Polynomial a)
list_to_polynomial p = case p of
   [x] -> x
   ((Const 0):xs) -> (list_to_polynomial xs)
   ((Const n):((Const m):xs)) -> list_to_monomial ( (Const (n+m)) : xs )
   (x:[y]) -> (add_monomial x y)
   (x:(y:xs)) 
    -> (if (eq_monomial x y) 
        then (PLUS (add_monomial x y) (list_to_polynomial xs))
        else (PLUS x (PLUS y (list_to_polynomial xs))))
 where
   eq_monomial x y = (not (gt_monomial x y)) && (not (gt_monomial y x))

----Normalizing polynomial
normalize :: (Num a) => (Polynomial a) -> (Polynomial a)
normalize p = list_to_polynomial (reverse (sort_monomial (flatten_polynomial (expand p))))

これだけだと使いづらいので、パーサ

import Text.ParserCombinators.Parsec
import Text.ParserCombinators.Parsec.Char
import Text.ParserCombinators.Parsec.Token
import Text.ParserCombinators.Parsec.Language
import Text.ParserCombinators.Parsec.Expr
import Monad
import Maybe
import Polynomial 


--字句解析
lexer  = makeTokenParser (haskellDef {reservedOpNames = ["+","*","-"]})
kwd   = reserved lexer
ident = identifier lexer

--構文解析
leftA s f = Infix (do { reservedOp lexer s; return f }) AssocLeft
table     = [[leftA "*" MULT] , [leftA "+" PLUS]]

polynomial_expr = buildExpressionParser table term
term = variable <|> constant <|> parens lexer polynomial_expr
constant    = liftM Const (natural lexer)
variable    = liftM Var ident

run :: Show a => Parser a -> String -> Maybe a
run p input = case (parse p "" input) of
              Left err -> Nothing
              Right x  -> Just x

pparse :: String -> Polynomial Integer
pparse str = fromJust (run polynomial_expr str)

Parsecはよく分からないし、なんでこれでいいのかもよく分かってない。とりあえず、動く

テスト。出力がS式なのは手抜き

*Main> (pparse "x*3 + y*(y+x*x+z)") + (pparse "x*x")
(+  (*  "x" (*  "x" "y")) (+  (*  "x" "x") (+  (*  3 "x") (+  (*  "y" "y") (*  "y" "z")))))

本当は、グレブナー基底の計算まで実装しようと思ってたけど、予想外に時間くったので、ここで終了。まあ、あと剰余計算書いたら、実質終わり

感想。なんというか、ひたすらパターンマッチだし、なんかS式書いてるのと、あんまり変わらないなー。別の言語というよりは、better Schemeな感じだ。

じゃあSchemeでこれを書いてたらどうなってたかというと、途中でめんどくさくなってたんじゃないかと。使用暦数年のScheme < 使用暦数日のHaskell。C/C++を使っていて、型があって助かると感じたことはあまりないけど、Haskellのは頭を整理する役に立つ。強い型付けも悪くないと思った。アセンブラをごりごりに書いたことがある人なら、Cの型でも便利と思ったりするんだろうか。実質上、Inductiveなデータ型というのは、S式とほとんど変わらないし、S式はInductiveなデータ構造の中で最も汎用的には違いないけど、Haskellみたいに型チェックしてくれる方が楽。

もっとも、数年前だったら挫折してたと思う。数年前は、Haskellについて書かれたサイトなんて数えるほどしかなかったし。今は、適当に検索すると、一杯ひっかかって、サンプルコードを適当に眺めれば大体なんとなく分かってくるので、実はそれが一番重要かもしれない

けど、最近、パターンマッチで済むことばっかやってるので、頭が悪くなるんじゃないかと心配。