如何在 Haskell 中计算正弦函数？

Question

这是我的问题：我需要一个 Haskell 函数，它使用关联的 Taylor serie 计算某个数字的 sine 的近似值...

我用 C++ 写了这个：

double msin(double number, int counter = 0, double sum = 0)
{
    // sin(x) = x - (x'3 / 3!) + (x'5 / 5!) - (x'7 / 7!) + (x'9 / 9!)
    if (counter <= 20)
    {
        if (counter % 2 == 0)
            sum += mpow(number, counter * 2 + 1) / mfak(counter * 2 + 1) ;
        else
            sum -= mpow(number, counter * 2 + 1) / mfak(counter * 2 + 1) ;

        counter++;
        sum =  msin(number, counter, sum);

        return sum;
    }

    return (sum* 180.0 / _PI);
}

现在我正在尝试在 Haskell 中执行此操作，但我不知道如何...现在我正在尝试这样的事情（它并没有真正起作用，但它正在进行中；））：

这行得通：

mfak number = if number < 2
                    then 1
                    else number *( mfak (number -1 )) 


mpow number potenca = if potenca == 0
                        then 0
                        else if potenca == 1
                        then 1
                        else (number * (mpow number (potenca-1)))

这不起作用：

msin :: Double -> Int -> Double -> Double                           
msin number counter sum = if counter <= 20
                                    then if counter `mod` 2==0
                                            then let sum = sum + (msin 1 (let counter = counter+1 in counter) sum) in sum 
                                            else let sum = sum + (msin 1 (let counter = counter+1 in counter) sum) in sum                                   
                                    else sum* 180.0 / 3.14

更新....无法编译：/“无法匹配预期类型Double' with actual type Int'”

msin :: Double -> Int -> Double -> Double                           
msin number counter sum = if counter <= 20
                                    then if counter `mod` 2==0
                                            then let sum' = sum + ((mpow number (counter*2+1))/(mfak counter*2+1)) in msin number (counter+1) sum'
                                            else let sum' = sum - ((mpow number (counter*2+1))/(mfak counter*2+1)) in msin number (counter+1) sum'                              
                                    else sum* 180.0 / 3.14

如您所见，最大的问题是如何在“sum”中添加一些东西，增加“counter”并使用这些新值再次递归......

P.S.我是 Haskell 的新手，所以请尽可能多地解释您的解决方案。我正在阅读一些教程，但是我找不到如何将某些表达式的结果保存到一个值中，然后在它之后继续其他代码......每次我尝试这样做时它都会返回我的值，并且我不想这样……

提前感谢您的帮助！

Answer 1

问题是像let stevec = stevec+1 in stevec 这样的表达式。 Haskell 不是命令式语言。这不会给stevec 加一。相反，它将stevec 定义为比自身多一的数字。不存在这样的数字，因此您将获得无限循环，或者，如果幸运的话，会崩溃。

而不是

stevec++;
vsota =  msin(stevilo, stevec, vsota);

你应该使用类似的东西

let stevec' = stevec + 1
in  msin stevilo stevec' vsota

或者只是

msin stevilo (stevec + 1) vsota

（这里还有一些我不明白的地方。您将需要mpow 和mfak。它们在哪里？）

Answer 2

我会稍微修改一下算法。首先我们可以定义阶乘逆列表：

factorialInv :: [Double]
factorialInv = scanl (/) 1 [1..]  -- 1/0! , 1/1! , 1/2! , 1/3! , ...

然后，我们遵循正弦系数：

sineCoefficients :: [Double]
sineCoefficients = 0 : 1 : 0 : -1 : sineCoefficients

然后，给定x，我们将上述两个列表乘以x 的幂，逐点：

powerSeries :: [Double]   -- ^ Coefficients
            -> Double     -- ^ Point x on which to compute the series
            -> [Double]   -- ^ Series terms
powerSeries cs x = zipWith3 (\a b c -> a * b * c) cs powers factorialInv
   where powers = iterate (*x) 1   -- 1 , x , x^2 , x^3 , ...

最后，我们取前 20 个术语并加以总结。

sine :: Double -> Double                   
sine = sum . take 20 . powerSeries sineCoefficients
 -- i.e., sine x = sum (take 20 (powerSeries sineCoefficients x))

Answer 3

As you can see the biggest problem is how to add something to "vsota",

在函数式语言中，您将在此处使用递归 - 变量 vstota 被实现为函数参数，在处理列表时从调用传递到调用。

例如，要对一个数字列表求和，我们可以这样写：

sum xs = go 0 xs
  where go total [] = total
        go total (x:xs) = go (total+x) xs

在命令式语言中，total 将是一个被更新的变量。这是一个函数参数，它被传递给go 的下一个递归调用。

在您的情况下，我将首先编写一个生成幂级数项的函数：

sinusTerms n x = ... -- the first n terms of x - (x'3 / 3!) + (x'5 / 5!) - (x'7 / 7!) ...

然后使用上面的sum函数：

sinus n x = sum (sinusTerms n x)

Answer 4

您也可以使用递归列表定义来获得[x, x^3, x^5 ...] 和[1, 1/3!, 1/5! ...] 无限序列。完成后，剩下的就是将它们的项目相乘并取和。

sinus count x = sum (take count $ zipWith (*) ifactorials xpowers)
    where xpowers     = x : map ((x*x)*) xpowers 
          ifactorials = 1 : zipWith (/) ifactorials [i*(i+1) | i <- [2, 4 .. ]]

另外，最好定义xpowers = iterate ((x*x)*) x，因为它似乎更具可读性。

Answer 5

我已尽力遵循您的约定。对于mfak 和mpow，您应该避免使用if，因为使用模式匹配来编写它们会更清晰：

mfak :: Int -> Int
mfak 0 = 1
mfak 1 = 1
mfak n = n * mfak (n - 1)

mpow :: Double -> Int -> Double
mpow _ 0 = 1
mpow x 1 = x
mpow x p = x * mpow x (p - 1)

在计算正弦之前，我们创建了一个系数列表[(sign, power, factorial)]：

x - (x^3 / 3!) + (x^5 / 5!) - (x^7 / 7!) + (x^9 / 9!)
→ [(1,1,1), (-1,3,6), (1,5,120), (-1,7,5040), (1,9,362880)]

列表是由列表推导式无限创建的。首先我们压缩列表[1,-1,1,-1,1,-1...] 和[1,3,5,7,9,11...]。这给了我们列表[(1,1), (-1,3), (1,5), (-1,7)...]。从这个列表中，我们创建了最终列表[(1,1,1), (-1,3,6), (1,5,120), (-1,7,5040)...]：

sinCoeff :: [(Double, Int, Double)]
sinCoeff = [ (fromIntegral s, i, fromIntegral $ mfak i)
           | (s, i) <- zip (cycle [1, -1]) [1,3..]]

（cycle 无限重复一个列表，[1,3..] 创建一个无限列表，从 1 开始，步长为 2）

最后，msin 函数接近定义。它还使用列表解析来实现其目标（请注意，我保留了* 180 / pi，尽管我不确定它是否应该在那里。Haskell 知道 pi）。

msin :: Int -> Double -> Double
msin n x = 180 * sum [ s * mpow x p / f | (s, p, f) <- take n sinCoeff] / pi

（take n sinCoeff 返回列表的第一个 n 元素）

您可以尝试使用以下代码：

main = do
    print $ take 10 sinCoeff
    print $ msin 5 0.5
    print $ msin 10 0.5

Answer 6

表达式的形式为 x*P(x²)。

为了获得最大效率，x² 中的多项式必须使用the Horner rule 计算，而不是单独计算 x² 的幂。

具有阶乘值的系数序列可以在 Haskell 中递归地表示，就像通常对斐波那契数列所做的那样。使用ghci 解释器作为我们的测试平台，我们有：

$ ghci
GHCi, version 8.8.4: https://www.haskell.org/ghc/  :? for help
 λ> 
 λ> 
 λ> nextCoeffs d c = c : (nextCoeffs (d+1) ((-c)/(fromIntegral $ (2*d+2)*(2*d+3))))
 λ> 
 λ> allCoeffs = nextCoeffs 0 1.0
 λ> 

其中 d 是系列内部的深度，c 是当前系数。

完整性检查：深度 3 处的系数必须是 7 的倒数！

 λ> 
 λ> 1.0 /(allCoeffs !! 3)
 -5040.0
 λ> 

Horner 规则可以通过 foldr1 :: (a -> a -> a) -> [a] -> a 库函数在 Haskell 中呈现。

按照 Haskell 的惯例，我冒昧地将术语 count 作为最左边的参数，因为它最有可能保持不变。这是出于柯里化（部分评估）的目的。

所以我们有：

 λ> :{
|λ>  msin count x = let { s = x*x ;  cs = take count allCoeffs ;
|λ>                       stepFn c acc = acc*s + c ; }
|λ>                 in  x * (foldr1 stepFn cs)
|λ> :}

健全性检查，包含 20 个术语：

 λ> 
 λ> pi
 3.141592653589793
 λ>
 λ> msin 20 (pi/6)
 0.49999999999999994
 λ> 
 λ> msin 20 (pi/2)
 1.0
 λ> 

旁注 1：最终乘以 180 / π 仅对 逆三角函数有用。

旁注 2：在实践中，为了获得相当快的收敛，应该使用正弦函数的周期性将输入变量 x 减少到 [-π,+π] 区间。