Scripts for paper: M. Ragni - "Mr.CAS - A Minimalistic (pure) Ruby CAS for Fast Prototyping and Code Generation"

00_README.md

Mr.CAS - A Minimalistic (pure) Ruby CAS for Fast Prototyping and Code Generation

Abstract

There are complete Computer Algebra System (CAS) systems on the marketm with complete solutions for analysis of analytical models. But exporting a model, for optimization or any other research activity, requires a lot of work, even with a good software.

This work presents a Ruby library that exposes minimalistic CAS capabilities - i.e. simplifications, substitutions, evaluations, etc. Library aims at rapid prototyping of numerical interfaces and code generation for different target languages, separating mathematical expression from exportation rules - e.g. models from numerical conditioning best practices.

The library is implemented in pure Ruby language and compatible with all Ruby interpreter flavours.

Gist Repository

This Gist repository contains a working version of the scripts presented in the different listings of the paper. Due to restriction on number of words, some obvious part were omitted. All the scripts complete are presented here, numbered as in the paper.

Expected outputs in console are highlighted with =>

Tasks

There is a Rakefile that runs the scripts. For installing dependencies:

License for Gist Repository

Copyright (c) 2016 Matteo Ragni

Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

This license covers all the file included in this repository

02_Listing_01.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

z  = CAS.vars 'z'        # creates a variable

f  = z ** 2 + 1          # define a symbolic expression
df = f.diff(z)           # derivative w.r.t. z
puts df                  # => (((z)^((2 - 1)) * 2 * 1) + 0)

g  = CAS.declare :g, f   # creates implicit expression
dg = g.diff(z)           # derivative w.r.t. z
puts dg                  # => ((((z)^((2 - 1)) * 2 * 1) + 0) * Dg[0](((z)^(2) + 1)))

exit 0

03_Listing_02.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

x, y = CAS::vars 'x', 'y'           # creates two variables

f = CAS.log( CAS.sin( y ) )         # symbolic expression
f.subs y => CAS.asin(CAS.exp(x))    # performs substitution
fs = f.simplify                     # simplifies expression
puts fs                             # => x

exit 0

04_Listing_03.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

x  = CAS.vars 'x'       # creates a variable
f  = x ** 2 + 1         # defines a symbolic expression
fc = f.call x => 2      # evaluates for x = 2
puts fc                 # => 5.0

exit 0

05_Listing_04.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

x, y = CAS.vars 'x', 'y'     # creates two variables
f    = CAS.declare :f, x     # declares implicit function :f
g    = CAS.declare :g, x, y  # declares implicit function :g

fg = f.greater_equal g       # :f is greater than :g
puts fg                      # => (f(x) ≥ g(x, y))

max_fg = CAS::max f, g
puts max_fg                  # => ((f(x) ≥ g(x, y)) ? f(x) : g(x, y))

exit 0

06_Listing_05.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

x = CAS::vars 'x'                    # creates a variable
f = CAS::log(CAS::sin(x))            # defines a symbolic function

proc = f.as_proc                     # compiles callable closure
r    = proc.call('x' => Math::PI/2)  # runs callable closure
puts r                               # => 0.0

exit 0

07_Listing_06.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

# Rules definition for Fortran Language
module CAS

  {
    # . . .
    CAS::Variable => Proc.new { "#{name}" },               # rule for Variable operation
    CAS::Sin      => Proc.new { "sin(#{x.to_fortran})" },  # rule for Sin operation
    # . . .
  }.each do |cls, prc|
    cls.send(:define_method, :to_fortran, &prc)            # Appends new to_fortran methods
  end

end

# Usage
x    = CAS.vars 'x'              # creates a variable
code = (CAS.sin(x)).to_fortran   # converts to target language
puts code                        # => sin(x)

exit 0

08_Example_01.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

%w|Mr.CAS
   Mr.CAS/c-opt|.each do |gem|
  raise LoadError, "Cannot load #{gem}" unless require gem
end

# Model
x, y = CAS.vars :x, :y                         # defines two variables
g    = CAS.declare :g, x                       # declares implicit function g(x)

f    = x ** y + g * CAS.log(CAS.sin(x ** y))   # f(x,y) = x^y + g(x) * log(sin(x^y))

# Code Generation
g.c_name = 'g_impl'                            # g(x) implementation token

clib = CAS::CLib.create "example" do
  include_local "g_impl"                       # includes local g_impl header
  implements_as "f_impl", f                    # token in interface of f(x,y)
end

puts clib.header
=begin
=>
// Header file for library: example.c

#ifndef example_H
#define example_H

// Standard Libraries
#include <math.h>

// Local Libraries
#include "g_impl"

// Definitions


// Functions
double f_impl(double x, double y);

#endif // example_H
=end

puts clib.source
=begin
=>
// Source file for library: example.c

#include "example.h"

double f_impl(double x, double y) {
  double __t_0 = pow(x, y);
  double __t_1 = g_impl(x);
  double __t_2 = sin(__t_0);
  double __t_3 = log(__t_2);
  double __t_4 = (__t_1 + __t_3);
  double __t_5 = (__t_0 + __t_4);

  return __t_5;
}

// end of example.c
// Source file for library: example.c

#include "example.h"

double f_impl(double x, double y) {
  double __t_0 = pow(x, y);
  double __t_1 = g_impl(x);
  double __t_2 = sin(__t_0);
  double __t_3 = log(__t_2);
  double __t_4 = (__t_1 + __t_3);
  double __t_5 = (__t_0 + __t_4);

  return __t_5;
}

// end of example.c
=end

exit 0

09_Example_02.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

%w|Mr.CAS
   Mr.CAS/c-opt|.each do |gem|
  raise LoadError, "Cannot load #{gem}" unless require gem
end

x = CAS::vars :x

# Parameters can be defined as implicit functions with no arguments
a = CAS.declare "PARAM_A"


g = (CAS.sqrt(x + a) - CAS.sqrt(x)) + CAS.sqrt(CAS::Pi + x)

# Particular Code Generation for difference between square roots.
# If the exported finds a difference between two square roots, it perform
# normalization. Please note that written in this way (for a sake of brevity)
# breaks optimized code generation flow for the argument of the square roots.
module CAS
  class Diff
    alias :__to_c_impl_old :__to_c_impl

    def __to_c_impl(v)
      if @x.is_a? CAS::Sqrt and @y.is_a? CAS::Sqrt
        "(#{@x.x.__to_c(v)} + #{@y.x.__to_c(v)}) / " +
        "( #{@x.__to_c(v)} + #{@y.__to_c(v)} )"
      else
        self.__to_c_impl_old(v)
      end
    end
  end
end

clib = CAS::CLib.create "g_impl" do
  define "PARAM_A()", 1.0     # defines the token PARAM_A()
                              # 1.0 is an arbitrary value for PARAM_A
  define "M_PI", Math::PI     # Defines the token M_PI, with value of
                              # Pi taken from interpreter
  implements_as "g_impl", g   # token in interface of f(x,y)
end

puts clib.header
=begin
=>
// Header file for library: g_impl.c

#ifndef g_impl_H
#define g_impl_H

// Standard Libraries
#include <math.h>

// Local Libraries


// Definitions
#define PARAM_A() 1.0
#define M_PI 3.141592653589793

// Functions
double g_impl(double x);

#endif // g_impl_H
// Header file for library: g_impl.c

#ifndef g_impl_H
#define g_impl_H

// Standard Libraries
#include <math.h>

// Local Libraries


// Definitions
#define PARAM_A() 1.0
#define M_PI 3.141592653589793

// Functions
double g_impl(double x);

#endif // g_impl_H
=end

puts clib.source
=begin
=>
// Source file for library: g_impl.c

#include "g_impl.h"

double g_impl(double x) {
  double __t_0 = PARAM_A();
  double __t_1 = (x + __t_0);
  double __t_2 = sqrt(__t_1);
  double __t_3 = sqrt(x);
  double __t_4 = (__t_1 + x) / ( __t_2 + __t_3 );
  double __t_5 = (M_PI + x);
  double __t_6 = sqrt(__t_5);
  double __t_7 = (__t_4 + __t_6);

  return __t_7;
}

// end of g_impl.c
// Source file for library: g_impl.c

#include "g_impl.h"

double g_impl(double x) {
  double __t_0 = PARAM_A();
  double __t_1 = (x + __t_0);
  double __t_2 = sqrt(__t_1);
  double __t_3 = sqrt(x);
  double __t_4 = (__t_1 + x) / ( __t_2 + __t_3 );
  double __t_5 = (M_PI + x);
  double __t_6 = sqrt(__t_5);
  double __t_7 = (__t_4 + __t_6);

  return __t_7;
}

// end of g_impl.c
=end

puts g
# => ((√((x + PARAM_A())) - √(x)) + √((π + x)))

exit 0

10_Example_03.py

#!/usr/bin/env python3

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

import sympy
import math

##
# Evaluates integral of symbolic function f
# starting from a, ending in b, with step size
# h = (b - a)/n
# Returns a float
def integrate(f, a, b, n):
    h = (b - a)/n
    x = sympy.symbols('x')
    func = sympy.lambdify((x), f)

    sums = (func(a) +
            func(b)) / 2.0

    for i in range(1, n):
        sums += func(a + i*h)

    return sums * h

##
# Evaluates order using a function f,
# a starting point a and endpoint b with
# mesh size n
# Returns a float
def order(f, a, b, n):
    x = sympy.symbols('x')

    f_ab = sympy.Subs(f, (x), (b)).n()-\
           sympy.Subs(f, (x), (a)).n()
    df   = f.diff(x)
    f_1n = integrate(df, a, b, n)
    f_2n = integrate(df, a, b, 2 * n)

    return math.log(
      (f_ab - f_1n) /
      (f_ab - f_2n),
    2)


x = sympy.symbols('x')
f = sympy.atan(x)

print(order(f, -1.0, 1.0, 100))
# => 1.9999999974244451

exit(0)

10_Example_03.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

##
# Evaluates integral of symbolic function f
# starting from a, ending in b, with step size
# h = (b - a)/n
# Returns a Float
def integrate(f, a, b, n)
  h = (b - a) / n

  func = f.as_proc

  sum = ((func.call 'x' => a) +
        (func.call 'x' => b)) / 2.0

  for i in (1...n)
    sum += (func.call 'x' => (a + i*h))
  end
  return sum * h
end

##
# Evaluates order using a function f,
# a starting point a and endpoint b with
# mesh size n
# Returns a Float
def order(f, a, b, n)
  x = CAS.vars 'x'

  f_ab = (f.call x => b) -
         (f.call x => a)
  df   = f.diff(x).simplify
  f_1n = integrate(df, a, b, n)
  f_2n = integrate(df, a, b, 2 * n)

  return Math.log(
    (f_ab - f_1n) /
    (f_ab - f_2n),
  2)
end

x = CAS.vars 'x'  # Defines a variable
f = CAS.arctan x  # Defines arctan(x)

puts(order f, -1.0, 1.0, 100)              # => 1.9999999974244451

exit 0

11_Example_04.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

%w|Mr.CAS
   Mr.CAS/c-opt|.each do |gem|
  raise LoadError, "Cannot load #{gem}" unless require gem
end

$x, $y, $h = CAS::vars :x, :y, :h
# Evaluates the factorial value of a
# CAS::Constant
def fact(n); (n < 2 ? 1 : n * fact(n - 1)); end

# Evaluates the vector of derivatives that will be used
# in the next function.
def coeff(f, n)
  df = [f]
  for _ in 2..n
    df << df[-1].diff($x).simplify + (df[-1].diff($y).simplify * df[0])
  end
  return df
end

# Generates the symbolic graph of the Taylor step
# that can be used by an integrator.
def taylor(f, n)
  df = coeff(f, n)
  y = $y
  for i in 0...df.size
    y = y + (($h ** (i + 1))/(fact(i + 1)) * df[i])
  end
  return y.simplify
end

# In this case the generator is tested against the
# function f = x * y
f = $x * $y
# And it is exported in a C library
clib = CAS::CLib.create "taylor" do
  implements_as "taylor_step", taylor(f, 4)
end
puts clib.header
=begin
=>
#ifndef taylor_H
#define taylor_H

// Standard Libraries
#include <math.h>

// Local Libraries


// Definitions


// Functions
double taylor_step(double y, double h, double x);

#endif // taylor_H
=end
puts clib.source
=begin
=>
// Source file for library: taylor.c

#include "taylor.h"

double taylor_step(double y, double h, double x) {
  double __t_0 = (x + y);
  double __t_1 = (h + __t_0);
  double __t_2 = pow(h, 2.0);
  double __t_3 = (__t_2) / (2.0 + M_EPSILON);
  double __t_4 = (x + __t_0);
  double __t_5 = (y + __t_4);
  double __t_6 = (__t_3 + __t_5);
  double __t_7 = pow(h, 3);
  double __t_8 = (__t_7) / (6 + M_EPSILON);
  double __t_9 = (y + x);
  double __t_10 = pow(x, 2.0);
  double __t_11 = (__t_10 + 1.0);
  double __t_12 = (__t_11 + __t_0);
  double __t_13 = (__t_0 + __t_9 + __t_12);
  double __t_14 = (__t_8 + __t_13);
  double __t_15 = pow(h, 4);
  double __t_16 = (__t_15) / (24 + M_EPSILON);
  double __t_17 = (2.0 + y);
  double __t_18 = (2.0 + x);
  double __t_19 = (__t_18 + __t_0);
  double __t_20 = (y + __t_11);
  double __t_21 = (__t_19 + __t_20);
  double __t_22 = (x + __t_11);
  double __t_23 = (__t_18 + __t_22);
  double __t_24 = (__t_23 + __t_0);
  double __t_25 = (__t_17 + __t_21 + __t_24);
  double __t_26 = (__t_16 + __t_25);
  double __t_27 = (y + __t_1 + __t_6 + __t_14 + __t_26);

  return __t_27;
}

// end of taylor.c
=end

12_Example_05.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

# Taylor expansion example

# Factorial
def fact(n); n < 2 ? 1 : n * fact(n - 1); end

# Taylor function. Takes a function as input,
# extracts the argument, and use it to perform
# the Taylor expansion, where f is a function
# and n is the order of the series
def Taylor(f, n)
  # We can use this method to intercept
  # the argument of the function and create
  # the central point for the Taylor series
  x = f.args()[0]
  x0 = CAS::vars (x.name.to_s + "0")
  # We create the derivatives
  df = [f]
  n.times do; df << df[-1].diff(x).simplify; end
  # and finally we accumulate all the argument
  # of the series, using factorial function
  tf = 0
  for i in 0..n
    tf += (df[i].subs(x => x0) / fact(i) * ((x - x0) ** i)).simplify
  end
  # The usage of simplifications allows us to
  # reduce drastically the size of the summation
  return tf.simplify
end

# Let's test it
z = CAS::vars "z"
f = CAS::sin(z)
tf = Taylor(f, 3)
puts tf
# => (sin(z0) +                               # order 0
#    (cos(z0) * (z - z0)) +                   # order 1
#    ((-sin(z0)) / (2) * ((z - z0))^(2)) +    # order 2
#    ((-cos(z0)) / (6) * ((z - z0))^(3)))     # order 3

13_Example_06.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

raise LoadError, 'Please install Mr.CAS first' unless require 'Mr.CAS'

# This module contains an exporter for Latex language.
# The exporter is not complete, but shows the capabilities
# and the semplicity on how a  language exporter can be written

module CAS
  {
    # Terminal nodes
    CAS::Constant              => Proc.new { "#{x}" },
    CAS::Variable              => Proc.new { "#{name}" },
    CAS::PI_CONSTANT           => Proc.new { "\\pi" },
    CAS::INFINITY_CONSTANT     => Proc.new { "\\infty" },
    CAS::NEG_INFINITY_CONSTANT => Proc.new { "-\\infty" },
    # Base functions
    CAS::Sum    => Proc.new { "\\left( #{x.map(&:latex).join(" + ")} \\right)" },
    CAS::Diff   => Proc.new { "\\left( #{x.latex} - #{y.latex} \\right)" },
    CAS::Prod   => Proc.new { "\\left( #{x.map(&:latex).join(" ")} \\right)"},
    CAS::Pow    => Proc.new { "{#{x.latex}}^{#{y.latex}}" },
    CAS::Div    => Proc.new { "\\dfrac{#{x.latex}}{#{y.latex}}" },
    CAS::Sqrt   => Proc.new { "\\sqrt{#{x.latex}}" },
    CAS::Invert => Proc.new { "-#{x.latex}" },
    CAS::Abs    => Proc.new { "\\left| #{x.latex} \\right|" },

    # Trigonometric functions
    CAS::Sin    => Proc.new { "\\sin \\left( #{x.latex} \\right)" },
    CAS::Asin   => Proc.new { "\\arcsin \\left( #{x.latex} \\right)" },
    CAS::Cos    => Proc.new { "\\cos \\left( #{x.latex} \\right)" },
    CAS::Acos   => Proc.new { "\\arccos \\left( #{x.latex} \\right)" },
    CAS::Tan    => Proc.new { "\\tan \\left( #{x.latex} \\right)" },
    CAS::Atan   => Proc.new { "\\arctan \\left( #{x.latex} \\right)" },

    # Trascendent functions
    CAS::Exp    => Proc.new { "e^#{x.latex}" },
    CAS::Ln     => Proc.new { "\\log \\left( #{x.latex} \\right)" },

    # Box Conditions
    CAS::BoxConditionOpen        => Proc.new { "#{lower.latex} < #{x.latex} < #{upper.latex}" },
    CAS::BoxConditionClosed      => Proc.new { "#{lower.latex} \\leq #{x.latex} \\leq #{upper.latex}" },
    CAS::BoxConditionUpperClosed => Proc.new { "#{lower.latex} < #{x.latex} \\leq #{upper.latex}" },
    CAS::BoxConditionLowerClosed => Proc.new { "#{lower.latex} \\leq #{x.latex} < #{upper.latex}" },

    # Conditions
    CAS::Equal        => Proc.new { "#{x.latex} = #{y.latex}" },
    CAS::Smaller      => Proc.new { "#{x.latex} < #{y.latex}" },
    CAS::Greater      => Proc.new { "#{x.latex} > #{y.latex}" },
    CAS::SmallerEqual => Proc.new { "#{x.latex} \\leq #{y.latex}" },
    CAS::GreaterEqual => Proc.new { "#{x.latex} \\geq #{y.latex}" },

    # Piecewise
    CAS::Piecewise => Proc.new {
      "\\left\\{ " +
      "  \\begin{array}{lr} " +
      "    #{x.latex} & #{condition.latex} \\\\" +
      "    #{y.latex}" +
      "  \\end{array}" +
      "\\right."
    },
    CAS::Max => Proc.new { "\\max \\left( #{x.latex}, \\, #{y.latex} \\right)" },
    CAS::Min => Proc.new { "\\min \\left( #{x.latex}, \\, #{y.latex} \\right)" },
    CAS::Function => Proc.new { "#{name}\\left( #{x.map(&:latex).join(", ")} \\right)" }


  }.each do |cls, blk|
    cls.send(:define_method, "latex", &blk)
  end
end

x, y = CAS::vars "x", "y"
f = CAS::declare "f", x, y
g = CAS::declare "g", x, y

h = CAS::Piecewise.new(
  f + g ** 2,
  g - CAS::sqrt(f),
  x.greater_equal(y))

puts h.latex

14_Example_07.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md

%w|Mr.CAS
   Mr.CAS/c-opt|.each do |gem|
  raise LoadError, "Cannot load #{gem}" unless require gem
end

# In this example a new Op (that implements the Sign function)
# will be implemented. A rule for the exportation in C is derived
# and also a version optimized with regularization that overloads
# the first one.

module CAS
  class Sign < CAS::Op
    # The fisrt function that is needed is the function
    # that is able to evaluate it
    def call(f)
      (@x.call(f) >= 0 ? 1 : -1)
    end

    # Now we formulate a derivative for the sign function
    # The derivative is quite questionable, and may be defined
    # using delta of dirac, or defined by considering
    # sign(x) = |x|/x for x != 0. We will go with this
    # second alternative.
    def diff(z)
      dx = @x.diff(z)
      (dx)/(abs(@x)) * (1 - (abs(dx) ** 2)/(@x ** 2))
    end

    # It is possible to give also a simplification
    # routine. For the sign simplification. The only simplification
    # we need is related to constant values. If we
    # have a constant, we can decide to eliminate the sign function
    def simplify
      @x = @x.simplify
      if @x.is_a? CAS::Constant
        return @x if @x.call({}) >= 0
        return -@x
      end
      return self
    end

    # Now we present a string representation for the
    # new op
    def to_s; "sign(#{@x})"; end

    # And a code representation that can ge interpreted by Ruby
    def to_code
      "(#{@x.to_code} >= 0 ? 1 : -1)"
    end
  end # end

  # It is possible to add a shortcut for calling rapidly the
  # op we just created
  def self.sign(x); CAS::Sign.new(x); end
end

# As additional step we can decide to append to our C exporter
# the sign function. Since it is only one Op to be added, we can go
# with a simpler solution with respect to the one presented in
# other examples. We also define the OPTIMIZED version, that requires
# three functions
module CAS
  class Sign
    # This function is the actual implementation
    # for the C code of the Op
    def __to_c_impl(v)
      "(#{@x.__to_c(v)} ? 1 : -1)"
    end

    # This function is responsible of the
    # creation of a new variable if does not exist.
    # it also creates the identification
    # The internal variable are always referred with __t_ID.
    def __to_c(v)
      sym = self.to_s.to_sym
      v[sym] = {
        def: "#{self.__to_c_impl(v)}",
        var: "__t_#{v.keys.size}",
        id: v.keys.size
      } unless v[sym]
      v[sym][:var]
    end

    # This is the user interface that is called by user.
    def to_c
      v = {}; self.__to_c(v); v
    end
  end
end

if ($0 == __FILE__)

  x, y = CAS::vars "x", "y"
  f = CAS.sign(x ** y)

  clib = CAS::CLib.create "sign_example" do
    implements_as "sign_example", f
  end

  puts clib.header
=begin
=>
// Header file for library: sign_example.c

#ifndef sign_example_H
#define sign_example_H

// Standard Libraries
#include <math.h>

// Local Libraries


// Definitions


// Functions
double sign_example(double x, double y);

#endif // sign_example_H
=end
  puts clib.source
=begin
=>
// Source file for library: sign_example.c

#include "sign_example.h"

double sign_example(double x, double y) {
  double __t_0 = pow(x, y);
  double __t_1 = (__t_0 ? 1 : -1);

  return __t_1;
}

// end of sign_example.c

=end
end

15_Example_08.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md


if ($0 == __FILE__)
  load "14_Example_07.rb"

  # Let's try to overload the sign op that we have defined in the
  # example 14_Example_07.rb for exporting a stabilized version of
  # the sign(x) => sin(arctan(a * x))

  x = CAS::vars "x"
  f = CAS.sign(x)

  # This is the version not stabilized
  clib_a = CAS::CLib.create "sign_unstable" do
    implements_as "sign", f
  end

  # Let's write an exportation rule stabilized
  module CAS; class Sign
    def __to_c_impl(v)
      "sin(atan(ALPHA * #{x.__to_c(v)}))"
    end
  end; end

  clib_b = CAS::CLib.create "sign_stable" do
    define "ALPHA", 1e4
    implements_as "sign", f
  end

  # Another way can be the substitution of
  # the sign function with a new symbolic version
  # but this CHANGES the graph of f
  a = CAS.const 1e4
  f.subs(CAS.sign(x) => CAS.sin(CAS.arctan(a * x)))

  clib_c = CAS::CLib.create "sign_symbolic" do
    implements_as "sign", f
  end

  puts clib_a.header, clib_a.source
=begin
=>
// Header file for library: sign_unstable.c

#ifndef sign_unstable_H
#define sign_unstable_H

// Standard Libraries
#include <math.h>

// Local Libraries


// Definitions


// Functions
double sign(double x);

#endif // sign_unstable_H

// Source file for library: sign_unstable.c

#include "sign_unstable.h"

double sign(double x) {
  double __t_0 = sin(atan(ALPHA * x));

  return __t_0;
}

// end of sign_unstable.c
=end
  puts clib_b.header, clib_b.source
=begin
=>
// Header file for library: sign_stable.c

#ifndef sign_stable_H
#define sign_stable_H

// Standard Libraries
#include <math.h>

// Local Libraries


// Definitions
#define ALPHA 10000.0

// Functions
double sign(double x);

#endif // sign_stable_H

// Source file for library: sign_stable.c

#include "sign_stable.h"

double sign(double x) {
  double __t_0 = sin(atan(ALPHA * x));

  return __t_0;
}

// end of sign_stable.c
=end
  puts clib_c.header, clib_c.source
=begin
=>
// Header file for library: sign_symbolic.c

#ifndef sign_symbolic_H
#define sign_symbolic_H

// Standard Libraries
#include <math.h>

// Local Libraries


// Definitions


// Functions
double sign(double x);

#endif // sign_symbolic_H

// Source file for library: sign_symbolic.c

#include "sign_symbolic.h"

double sign(double x) {
  double __t_0 = sin(atan(ALPHA * x));

  return __t_0;
}

// end of sign_symbolic.c
=end
end

16_Example_09.rb

#!/usr/bin/env ruby

# License MIT - 2016 Matteo Ragni :: refer to 0_README.md


%w|Mr.CAS
   Mr.CAS/auto-diff|.each do |gem|
  raise LoadError, "Cannot load #{gem}" unless require gem
end

def NewtonRaphson(f, x, x0, config = {})
  # configurations for the algorithm
  cfg = {maxsteps: 40, tolerance: 1e-12}
  # Configuration overrides
  config.keys.each { |k| cfg[k] = config[k] }
  # Starting Point
  xv = x0

  0.upto(cfg[:maxsteps]) do
    fv = f.auto_diff({x => xv})
    # We return the result if we reached
    # the zero of the function
    return xv if fv.x.abs < cfg[:tolerance]
    xv = xv - fv.x/fv.y
  end
  # If we cannot reach a minimum we
  # raise an argument error
  raise RuntimeError, "Not Found"
end


x, y = CAS::vars "x", "y"
a, b = CAS::const 1.0, 100.0

f = ((a - x) ** 2) + b * ((y - (x ** 2)) ** 2)
f.subs({y => (x ** 2)}).simplify

x0 = NewtonRaphson(f, x, 0.0)
puts "#{f} = #{f.call({x => x0})} =~ 0 (for x = #{x0})"

Rakefile

desc "Install ruby dependencies (Mr.CAS.) for current user"
task :dep_ruby do
  `gem install Mr.CAS`
end

desc "Install ruby dependencies (Mr.CAS.) system wide"
task :dep_ruby_root do
  `sudo gem install Mr.CAS`
end

desc "Install python dependencies (sympy) for current user"
task :dep_py do
  `pip3 install sympy`
end

desc "Install python dependencies (sympy) system wide"
task :dep_py do
  `pip3 install sympy`
end

require 'pry-byebug'

task_list = []
Dir.glob("*.rb").sort.each do |script|
  name = script.gsub(/\d+_(List|Ex)\w+_(\d+)\.rb/, '\1\2')

  desc "Run Ruby script #{script}"
  task name.to_sym do
    puts "== Running Ruby task :: #{name} :: #{script} =="
    puts `ruby #{script}`
    puts
  end

  task_list << name.to_sym
end

desc "Run Python3 script 10_Example_03.py"
task :Ex03py do
  puts "== Running Python task :: Example03py :: 10_Example_03.py =="
  puts `python3 10_Example_03.py`
  puts
end

task_list << :Ex03py

desc "Runs all scripts in paper's order of appearance"
task :default => task_list