(**
How to formally differentiate a function
This small program shows how differentiation can be seen as a purely formal and automatic procedure.
Program in OCaml v4+ by Lilian Besson (ENS Cachan) for MA 101 at Mahindra École Centrale, (C) 2014.
GPLv3 Licenced.
@date Monday 14 September 2014
@author Lilian BESSON
*)
print_endline "# How to (simply) differentiate a function thanks to OCaml formal programming.\n";;
(** Define the data type.
*)
(** A function is an expression depending on some variables.
We choose to work only with [sin], [cos], [exp] and [ln] as basic functions. *)
type funct = Id of string | Lambda of float
| Sin of funct | Cos of funct | Exp of funct | Ln of funct
| Prod of funct * funct | Inv of funct | Sum of funct * funct
| Power of funct * int
| Comp of funct * funct
| Unknown of string
;;
(** Examples.
*)
let x = Id("x");;
let fsin = Sin(x)
and fcos = Cos(x)
and fexp = Exp(x)
and fln = Ln(x);;
let y = Id("y");;
let fpi = Lambda(3.1415)
and fpoly1 = Sum(Lambda(1.), Power(y, 7))
and fcossin = Comp((Cos(y)), (Sin(y)))
and fexpln = Comp((Exp(y)), (Ln(y)))
and finv = Inv(y)
and fsqsq = Power(y, 4)
and fh = Unknown("h");;
(** Simplify as much as possible.
*)
(** Expression returned if the function is not valid (some basic check, not sophiscated ones). *)
exception Not_defined;;
(** Recursively put in canonic form. *)
let rec canonic = function
| Sin(Lambda(0.)) -> Lambda(0.)
| Cos(Lambda(0.)) -> Lambda(1.)
| Ln(Lambda(x)) when x<=0. -> raise Not_defined
| Ln(Lambda(1.)) -> Lambda(0.)
| Prod(Lambda(0.), g) -> Lambda(0.)
| Prod(Lambda(1.), g) -> (canonic g)
| Prod(Lambda(-1.), Prod(Lambda(-1.), g)) -> g
| Prod(g, Lambda(x)) -> canonic (Prod(Lambda(x), (canonic g)))
| Prod(Lambda(l), g) -> Prod(Lambda(l), (canonic g))
| Prod(f,g) -> Prod((canonic f), (canonic g))
| Inv(Inv(f)) -> f
| Inv(Lambda(0.)) -> raise Not_defined
| Inv(Lambda(x)) -> Lambda(1. /. x)
| Inv(f) -> Inv((canonic f))
| Power(f,0) -> Lambda(1.)
| Power(f,1) -> canonic f
| Power(f,-1) -> canonic (Inv(f))
| Power(f,n) when n<0 -> canonic (Inv(Power(f, -n)))
| Power(f,n) -> Power((canonic f), n)
| Sum(f, Prod(Lambda(-1.), g) ) when f=g -> Lambda(0.)
| Sum(f,g) -> Sum((canonic f), (canonic g))
| Comp(Exp(Id(_)),Ln(f)) -> f
| Comp(Ln(f),Exp(Id(_))) -> f
| Comp(f,g) -> Comp((canonic f), (canonic g))
| u -> u
;;
(** Pretty printing a function.
*)
let print = Format.printf;;
let sprint = Format.sprintf;;
let fprint = Printf.fprintf;;
(** To a string without the variable name. *)
let rec string_of_funct_novar = function
| Id(x) -> ""
| Lambda(l) -> sprint "%g" l
| Sin(Id(_)) -> "sin"
| Cos(Id(_)) -> "cos"
| Exp(Id(_)) -> "exp"
| Ln(Id(_)) -> "ln"
| Sin(f) -> sprint "sin(%s)" (string_of_funct_novar f)
| Cos(f) -> sprint "cos(%s)" (string_of_funct_novar f)
| Exp(f) -> sprint "exp(%s)" (string_of_funct_novar f)
| Ln(f) -> sprint "ln(%s)" (string_of_funct_novar f)
| Prod(Lambda(0.), g) -> "0"
| Prod(Lambda(1.), g) -> (string_of_funct_novar g)
| Prod(Lambda(-1.), g) -> sprint "-(%s)" (string_of_funct_novar g)
| Prod(Lambda(l), g) -> sprint "%g(%s)" l (string_of_funct_novar g)
| Prod(f,g) -> sprint "(%s).(%s)" (string_of_funct_novar f) (string_of_funct_novar g)
| Inv(f) -> sprint "(1/(%s))" (string_of_funct_novar f)
| Sum(f,g) -> sprint "(%s) + (%s)" (string_of_funct_novar f) (string_of_funct_novar g)
| Power(f,0) -> "1"
| Power(f,1) -> (string_of_funct_novar f)
| Power(f,-1) -> sprint "(1/(%s))" (string_of_funct_novar f)
| Power(f,n) when n<0 -> sprint "(%s)^(%i)" (string_of_funct_novar f) n
| Power(f,n) -> sprint "(%s)^%i" (string_of_funct_novar f) n
| Comp(f,g) -> sprint "%s(%s)" (string_of_funct_novar f) (string_of_funct_novar g)
| Unknown(s) -> s
;;
(** To a string with the variable name. *)
let rec string_of_funct = function
| Id(x) -> x
| Lambda(l) -> sprint "%g" l
| Sin(f) -> sprint "sin(%s)" (string_of_funct f)
| Cos(f) -> sprint "cos(%s)" (string_of_funct f)
| Exp(f) -> sprint "exp(%s)" (string_of_funct f)
| Ln(f) -> sprint "ln(%s)" (string_of_funct f)
| Prod(Lambda(0.), g) -> "0"
| Prod(Lambda(1.), g) -> (string_of_funct g)
| Prod(Lambda(-1.), g) -> sprint "-(%s)" (string_of_funct g)
| Prod(Lambda(l), g) -> sprint "%g(%s)" l (string_of_funct g)
| Prod(f,g) -> sprint "(%s).(%s)" (string_of_funct f) (string_of_funct g)
| Inv(Inv(f)) -> (string_of_funct f)
| Inv(Id(x)) -> sprint "1/%s" x
| Inv(f) -> sprint "(1/(%s))" (string_of_funct f)
| Sum(Lambda(l),g) -> sprint "%g + (%s)" l (string_of_funct g)
| Sum(g,Lambda(l)) -> sprint "(%s) + %g" (string_of_funct g) l
| Sum(f,g) -> sprint "(%s) + (%s)" (string_of_funct f) (string_of_funct g)
| Power(f,0) -> "1"
| Power(f,1) -> (string_of_funct f)
| Power(f,-1) -> sprint "(1/(%s))" (string_of_funct f)
| Power(f,n) when n<0 -> sprint "(%s)^(%i)" (string_of_funct f) n
| Power(f,n) -> sprint "(%s)^%i" (string_of_funct f) n
| Comp(f,g) -> sprint "%s(%s)" (string_of_funct_novar f) (string_of_funct g)
| Unknown(s) -> s
;;
(** Printing. *)
let print_funct_novar f = print_endline (string_of_funct_novar (canonic f));;
(** Printing. *)
let print_funct f = print_endline (string_of_funct (canonic f));;
let print_funct_fancy name x f = print " %s : %s --> %s\n" name x (string_of_funct (canonic f));;
(** Examples.
*)
print "## First five functions :\n";;
print_funct_fancy "Identity" "x" x;;
print_funct_fancy "Sine" "x" fsin;;
print_funct_fancy "Cosine" "x" fcos;;
print_funct_fancy "Exponential" "x" fexp;;
print_funct_fancy "Logarithm" "x" fln;;
print "\n\n## Some other examples of functions :\n";;
print_funct_fancy "Identity" "y" y;;
print_funct_fancy "Constant" "y" fpi;;
print_funct_fancy "Poly1" "y" fpoly1;;
print_funct_fancy "cos o sin" "y" fcossin;;
print_funct_fancy "exp o ln" "y" fexpln;;
print_funct_fancy "Inverse" "y" finv;;
print_funct_fancy "Power 4" "y" fsqsq;;
print_funct_fancy "Unknown h" "y" fh;;
(** Differentiate with formal rules.
*)
let rec differentiate_aux = function
| Lambda(l) -> Lambda(0.)
| Id(x) -> Lambda(1.)
| Sin(f) -> Prod( (differentiate_aux f), Cos(f) )
| Cos(f) -> Prod( (differentiate_aux f), ( Prod(Lambda(-1.), (Sin(f)))) )
| Exp(f) -> Prod( (differentiate_aux f), Exp(f) )
| Ln(f) -> Prod( (differentiate_aux f), Inv(f) )
| Prod(f,g) -> Sum( (Prod((differentiate_aux f), g)), (Prod(f, (differentiate_aux g))) )
| Inv(f) -> Prod( ( Prod(Lambda(-1.), (differentiate_aux f)) ), (Power(f, -2)) )
| Sum(f,g) -> Sum((differentiate_aux f), (differentiate_aux g))
| Power(f,0) -> Lambda(0.)
| Power(f,1) -> differentiate_aux f
| Power(f,n) -> Prod( (Prod( (Lambda(float_of_int n)), (differentiate_aux f) )), (Power( f, n-1 )) )
| Comp(f,g) -> Prod( (differentiate_aux g), (Comp( (differentiate_aux f), g )) )
| Unknown(s) -> Unknown( (sprint "%s'" s) )
;;
(** Final value, when we put in simplify [f] as much as possible before differentiating. *)
let differentiate f = (canonic (differentiate_aux (canonic f)));;
(** Examples.
*)
print "\n\n## Differential of the first five functions :\n";;
print_funct_fancy "Differential of Identity" "x" (differentiate x);;
print_funct_fancy "Differential of Sine" "x" (differentiate fsin);;
print_funct_fancy "Differential of Cosine" "x" (differentiate fcos);;
print_funct_fancy "Differential of Exponential" "x" (differentiate fexp);;
print_funct_fancy "Differential of Logarithm" "x" (differentiate fln);;
print "\n\n## Differential of some other examples of functions :\n";;
print_funct_fancy "Differential of Identity" "y" (differentiate y);;
print_funct_fancy "Differential of Constant" "y" (differentiate fpi);;
print_funct_fancy "Differential of Poly1" "y" (differentiate fpoly1);;
print_funct_fancy "Differential of cos o sin" "y" (differentiate fcossin);;
print_funct_fancy "Differential of exp o ln" "y" (differentiate fexpln);;
print_funct_fancy "Differential of Inverse" "y" (differentiate finv);;
print_funct_fancy "Differential of Power 4" "y" (differentiate fsqsq);;
print_funct_fancy "Differential of Unknown h" "y" (differentiate fh);;
(** Fin *)
print "\n\nEnd of the tests, program in OCaml v4+ by Lilian Besson (ENS Cachan), for MA 101 at Mahindra École Centrale, (C) 2014.";;
print "\n~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n";;