(module math.flonum.log1pmx (fplog1pmx)
  (import scheme
          chicken.type
          math.flonum.functions
          math.flonum.polyfun
          math.flonum.log
          math.polynomial.chebyshev)

  (define fplog1pmx-neg-0.25-0.5
    (inline-chebyshev-fppoly-fun
     -0.5 -0.25
     (-0.2103141053572070997362823489693816811996
      0.07704102886728760721086370117643443131217
      -0.01020514433643803605431850588217215237178
      6.872859538243712918354721556695090979546e-4
      -5.207248546376666160678822643500842976255e-5
      4.208311415510517753022518076392447733423e-6
      -3.542714867432068706621663292676355004748e-7
      3.067601814854621065962689625154653762101e-8
      -2.711543742374190275109136442040373034036e-9
      2.434858166790830434687283890247233106842e-10
      -2.21373562123951709989617540939460847654e-11
      2.033024647979063139529426446529956536822e-12
      -1.882624294787569740261613339485357785274e-13
      1.755541612919310465016807409127750031632e-14
      -1.646781645171118497659636869130525525145e-15
      1.552680302240995074747799623556383281355e-16
      -1.470371377726707214575139072206825846451e-17
      1.385352010971336691234958236784925394612e-18)))

  (: fllog1pmx-taylor-0 (float -> float))
  (define (fplog1pmx-taylor-0 x)
    (fp* (fp* x x)
         ((make-fppolyfun
           (#i-1/2 #i1/3 #i-1/4 #i1/5 #i-1/6 #i1/7 #i-1/8 #i1/9 #i-1/10 #i1/11
                   #i-1/12 #i1/13 #i-1/14 #i1/15 #i-1/16 #i1/17 #i-1/18 #i1/19 #i-1/20 #i1/21
                   #i-1/22 #i1/23 #i-1/24 #i1/25 #i-1/26))
          x)))

  (define fplog1pmx-0.25-0.5
    (inline-chebyshev-fppoly-fun
     0.25 0.5
     (-0.1172376347903564395048972641928223555326
      -0.03390230020664453827879131203208535815299
      -0.002074697726910078933295567647061059776673
      6.300006356267044793209616759065448160092e-5
      -2.15218532902292423180294171512095168131e-6
      7.842365320119772047407362890838076306055e-8
      -2.976756006675387489388663469498613703626e-9
      1.16218125023220889846755872728148163628e-10
      -4.631901690461512483973214435928711719421e-12
      1.875358176311070533829034898220631704307e-13
      -7.687836726698030051212239319592396588335e-15
      3.183382900679105044910163409271165546662e-16
      -1.329157157589045915049329706194501559878e-17
      5.578423096000352584032742755315719988568e-19)))

  (: fplog1pmx-big (float -> float))
  (define (fplog1pmx-big x)
    (cond [(fp< (fpabs x) 1.0)
           (let ((y (fp+ x 1.0)))
             (fp- (fp- (fplog y) x) (fp/ (fp- (fp- y 1.0) x) y)))]
          [else
           (fp- (fplog1p x) x)]))

  (: fplog1pmx (float -> float))
  (define (fplog1pmx x)
    (cond [(fp< x 0.0)
           (cond [(fp> x -0.25)  (fplog1pmx-taylor-0 x)]
                 [(fp> x -0.50)  (fplog1pmx-neg-0.25-0.5 x)]
                 [(fp> x -1.00)  (fplog1pmx-big x)]
                 [(fp= x -1.00)  -inf.0]
                 [else  +nan.0])]
          [(fp> x 0.0)
           (cond [(fp< x 0.25)  (fplog1pmx-taylor-0 x)]
                 [(fp< x 0.50)  (fplog1pmx-0.25-0.5 x)]
                 [else  (fplog1pmx-big x)])]
          [else  0.0]))
  )