search

mathnet / mathnet-symbolics

Math.NET Symbolics
http://symbolics.mathdotnet.com
MIT License
341 stars 66 forks source link

algebraically exact values of sine function #44

Open diluculo opened 6 years ago

diluculo commented 6 years ago

This Symbolics says that "Expressions always appear in a simplified form according to a set of rules."

If so, how about adding some known rules to the trigonometric functions? For some m/n*pi (with integer m and n), trigonometric functions can be simply expressed as algebraic numbers. For example, sin(pi/3) returns sqrt(3)/2. However, we can't handle all of cases, so, if we limit to just a few cases, we can easily apply this rule with lookup tables.

Rules can be added: (1) Argument that are rational multiple of pi can be expressed in algebraic numbers: sin(11/12*pi) = (sqrt(3) - 1)/(2*sqrt(2)), ... (2) Hence, exponential function can be rewritten as cos and sin terms: exp(m/n*pi*j) = cos(m/n*pi) + j*sin(m/n*pi), where j is imaginary unit (3) Argument that are multiple of j can be expressed in terms of hyperbolic functions: sin(z*j) = j*sinh(z), ...

diluculo commented 6 years ago

There are many rules we can apply to elementary functions such as exp(x), sin(x) and sinh(x). Here, I would like to ask which of the following is appropriate for the purpose of this library. I'd like to send a PR for selected rules.

(1) auto-simplification of argument: sin(x + y - (x + y)) = sin(0) = 0 (2) undefined and infinities : sin(infinity) = undefined (3) zero, one, pi, and j: sin(pi) = 0, sin(j) = j*sinh(1) (4) known rational multiples of pi: sin(1/4*pi) = 1/sqrt(2), ... (5) multiples of j: sin(j*x) = j*sinh(x) (6) remove pi term: sin(2*pi + x) = sin(x), sin(pi + x) = -sin(x), sin(1/2*pi + x) = cos(x), sin(3/2*pi + x) = -cos(x) (7) negative argument: sin(-x) = -sin(x) (8) negative first term: sin(-x + y) = -sin(x - y) (9) forward-inverse identities: sin(asin(x)) = x and more....

(2), (3) and (7) are already partially implemented.

cdrnet commented 6 years ago

Regarding (1): in general, auto-simplification should not expand expressions. However, in this case it might be a bug that the inner expression is not auto-simplified to zero already (need to rethink this). So, once fixed, this should work out of the box.

diluculo commented 6 years ago

@cdrnet, I think so. Now it is clear -(a + b) is not auto-simplified to -a - b. In that respect, returns in the below table are very acceptable, except (d). For consistency, how about to allow only minus or -1 to be distributed inside the parenthesis. i.e. -(a + b) = -a - b. I think the last column may be the correct behavior.

No Input Return Expected
(a) a + b - (a + b) a + b - (a + b) 0
(b) a + b + (-a - b) 0 0
(c) a + b - 2*(a + b) a + b - 2*(a + b) a + b - 2*(a + b)
(d) -(a + b) + 2*(a + b) a + b -a - b + 2*(a + b)
(e) -a - b + 2*(a + b) -a - b + 2*(a + b) -a - b + 2*(a + b)
(f) 2(a + b) - 2(a + b) 0 0
cdrnet commented 6 years ago

Yes, we'd only distribute -1, otherwise the product is not just a technical artifact of our internal representation (-x is (-1)*x internally). d) would no longer auto-simplify to a+b, but that may be ok.

I agree, I'd also expect the auto-simplified form of the last column.

diluculo commented 6 years ago

I can get the last column by simply adding a rule(| Sum ax as x' ...) to the valueMul function as shown below. 

and multiply x y =
      ...
      let rec valueMul (v:Value) x =
            if Value.isZero v then zero else
            match x with
            ...
            | Product ax -> if Value.isOne v then x else Product (Values.unpack v::ax)
            | Sum ax as x'-> if Value.isOne v then x' // intentionally distribute -1, e.g. -(a + b) = -a - b
                             else if Value.isMinusOne v then ax |> List.map (function xi -> valueMul v xi) |> Sum
                             else Product [Values.unpack v; x']
            | x -> if Value.isOne v then x else Product [Values.unpack v; x]