The low-level interface

Mikkel Meyer Andersen and Søren Højsgaard

2023-01-16

library(Ryacas)

The low-level interface consists of these two main functions:

Note, that the yacas command x is a string and must often be built op using paste()/paste0(). Examples of this will be shown in multiple examples below.

A short summary of often-used yacas commands are found in the section “yacas reference” in the “Getting started” vignette. A short summary of Ryacas’s low-level functions are also found in the section “Ryacas low-level reference” at the end of this document.

Note that the yacas documentation is a very useful resource.

Example 1: Simple algebra

First, consider this polynomial:

eq <- "x^2 + 4 + 2*x + 2*x"

Now, perform yacas operations, and get result as string/character:

yac_str(eq) # No task was given to yacas, so we simply get the same returned
## [1] "x^2+2*x+2*x+4"
yac_str(paste0("Simplify(", eq, ")"))
## [1] "x^2+4*x+4"
yac_str(paste0("Factor(", eq, ")"))
## [1] "(x+2)^2"
yac_str(paste0("TeXForm(Factor(", eq, "))"))
## [1] "\\left( x + 2\\right)  ^{2}"

\[ \left( x + 2\right) ^{2} \]

Also see yacas documentation on Simplify(), Factor() and TeXForm().

Instead of the pattern paste0("Simplify(", eq, ")") etc., there exists a helper function y_fn() that does this:

y_fn(eq, "Simplify")
## [1] "Simplify(x^2 + 4 + 2*x + 2*x)"
yac_str(y_fn(eq, "Simplify"))
## [1] "x^2+4*x+4"
yac_str(y_fn(eq, "Factor"))
## [1] "(x+2)^2"
yac_str(y_fn(y_fn(eq, "Factor"), "TeXForm"))
## [1] "\\left( x + 2\\right)  ^{2}"

As you see, there are a lot of nested function calls. That can be avoided by using magrittr’s pipe %>% (automatically available with Ryacas) together with the helper function y_fn():

eq %>% y_fn("Simplify")
## [1] "Simplify(x^2 + 4 + 2*x + 2*x)"
eq %>% y_fn("Simplify") %>% yac_str()
## [1] "x^2+4*x+4"
eq %>% y_fn("Factor") %>% yac_str()
## [1] "(x+2)^2"
eq %>% y_fn("Factor") %>% y_fn("TeXForm") %>% yac_str()
## [1] "\\left( x + 2\\right)  ^{2}"

Below, we will stick to the standard way of calling the functions, and not using the pipe, but now it has been demonstrated if the user prefers that way.

We will not use the pipe operator below, but just demonstrate its usage.

Now, again perform yacas operations, but get result as an R expression, e.g. for continued computations:

eq
## [1] "x^2 + 4 + 2*x + 2*x"
eq %>% yac_expr() # Alternative to "yac_expr(eq)"
## expression(x^2 + 2 * x + 2 * x + 4)
cmd <- eq %>% y_fn("Factor")
cmd
## [1] "Factor(x^2 + 4 + 2*x + 2*x)"
e <- yac_expr(cmd)
e
## expression((x + 2)^2)
eval(e, list(x = 2))
## [1] 16

Example 2: Linear algebra

To work with matrices and vectors, you need to realise that yacas and R has different ways of representing these objects. yacas represents vectors as a list, and a matrix as a list of lists (each list is a row).

Simple example

You can work with these directly. Here illutrated with vectors:

cmd <- "2*{x, x^2, x^3}"
cmd %>% yac_str()
## [1] "{2*x,2*x^2,2*x^3}"
e <- cmd %>% yac_expr()
e
## expression(c(2 * x, 2 * x^2, 2 * x^3))
eval(e, list(x = 1.5))
## [1] 3.00 4.50 6.75

And then illutrated with matrices. First with purely numeric contents:

cmd <- "{{1, 2}, {3, 4}}"
yac_str(cmd)
## [1] "{{1,2},{3,4}}"
y_print(cmd) # Convenience function for prettier print
## {{1, 2},
##  {3, 4}}
e <- cmd %>% yac_expr()
e
## expression(rbind(c(1, 2), c(3, 4)))
eval(e)
##      [,1] [,2]
## [1,]    1    2
## [2,]    3    4

Also in \(\LaTeX\) (yacas documentation on TeXForm()):

cmd %>% y_fn("TeXForm") %>% yac_str()
## [1] "\\left( \\begin{array}{cc} 1 & 2 \\\\ 3 & 4 \\end{array} \\right) "

But it also works with symbolic contents:

cmd1 <- paste0("a * ", cmd, "")
cmd2 <- cmd1 %>% y_fn("Inverse")
cmd2
## [1] "Inverse(a * {{1, 2}, {3, 4}})"
cmd2 %>% yac_str()
## [1] "{{1/a+(6*a^2)/(a^2*(4*a-(6*a^2)/a)),(-(2*a)/a)/(4*a-(6*a^2)/a)},{(-(3*a)/a)/(4*a-(6*a^2)/a),1/(4*a-(6*a^2)/a)}}"
cmd2 %>% y_fn("TeXForm") %>% yac_str()
## [1] "\\left( \\begin{array}{cc} \\frac{1}{a}  + \\frac{6 a ^{2}}{a ^{2} \\left( 4 a - \\frac{6 a ^{2}}{a} \\right) }  & \\frac{ - \\frac{2 a}{a} }{4 a - \\frac{6 a ^{2}}{a} }  \\\\ \\frac{ - \\frac{3 a}{a} }{4 a - \\frac{6 a ^{2}}{a} }  & \\frac{1}{4 a - \\frac{6 a ^{2}}{a} }  \\end{array} \\right) "

\[ \left( \begin{array}{cc} \frac{1}{a} + \frac{6 a ^{2}}{a ^{2} \left( 4 a - \frac{6 a ^{2}}{a} \right) } & \frac{ - \frac{2 a}{a} }{4 a - \frac{6 a ^{2}}{a} } \\ \frac{ - \frac{3 a}{a} }{4 a - \frac{6 a ^{2}}{a} } & \frac{1}{4 a - \frac{6 a ^{2}}{a} } \end{array} \right) \]

paste0(cmd2, "*", cmd1) %>% 
  y_fn("Simplify") %>% 
  yac_str()
## [1] "{{1,0},{0,1}}"
e2 <- cmd2 %>% yac_expr()
eval(e2, list(a = 2.2))
##            [,1]       [,2]
## [1,] -0.9090909  0.4545455
## [2,]  0.6818182 -0.2272727

Using R’s character matrices

The above is fine when writing yacas vectors and matrices by hand. But often one would want to exploit R’s convenient functions to work with matrices.

The central idea to make this possible is to work with R character matrices. We provide two helper functions to go back and forth between R and yacas:

Below, we illustrate the usage of both functions.

First, we create a character matrix using R:

Achr <- matrix(0, nrow = 3, ncol = 3)
diag(Achr) <- 1
Achr[2, 3] <- "a"
Achr[1, 3] <- "a"
Achr
##      [,1] [,2] [,3]
## [1,] "1"  "0"  "a" 
## [2,] "0"  "1"  "a" 
## [3,] "0"  "0"  "1"

Note how this is a character matrix. If we want to find it’s inverse symbolically using yacas, it must first be represented as a yacas matrix:

Ayac <- Achr %>% as_y()
Ayac
## [1] "{{1, 0, a}, {0, 1, a}, {0, 0, 1}}"

Now, we can find the inverse:

cmd <- Ayac %>% y_fn("Inverse")
cmd %>% yac_str()
## [1] "{{1,0,-a},{0,1,-a},{0,0,1}}"

A nicer representation can be obtained in (at least) four ways:

way1 <- cmd %>% yac_str()
way1 %>% y_print()
## {{ 1,  0, -a},
##  { 0,  1, -a},
##  { 0,  0,  1}}
way2 <- cmd %>% y_fn("TeXForm") %>% yac_str()
way2
## [1] "\\left( \\begin{array}{ccc} 1 & 0 &  - a \\\\ 0 & 1 &  - a \\\\ 0 & 0 & 1 \\end{array} \\right) "
way3 <- cmd %>% y_fn("PrettyForm") %>% yac_str()
way3 %>% cat() # Result of PrettyForm() must be printed
## 
## /                         \
## | ( 1 ) ( 0 ) ( -( a ) )  |
## |                         |
## | ( 0 ) ( 1 ) ( -( a ) )  |
## |                         |
## | ( 0 ) ( 0 ) ( 1 )       |
## \                         /
way4 <- cmd %>% yac_str()
way4
## [1] "{{1,0,-a},{0,1,-a},{0,0,1}}"
way4 %>% as_r()
##      [,1] [,2] [,3]
## [1,] "1"  "0"  "-a"
## [2,] "0"  "1"  "-a"
## [3,] "0"  "0"  "1"
way4 %>% as_r() %>% print(quote = FALSE)
##      [,1] [,2] [,3]
## [1,] 1    0    -a  
## [2,] 0    1    -a  
## [3,] 0    0    1

Say we want to subset it to only consider a submatrix. To do that, we can use R’s facilities:

A_inv_yac <- way4 %>% as_r()
Bchr <- A_inv_yac[2:3, 2:3]
Bchr
##      [,1] [,2]
## [1,] "1"  "-a"
## [2,] "0"  "1"
Bchr %>% as_y()
## [1] "{{1, -a}, {0, 1}}"
Bchr %>% as_y() %>% 
    y_fn("Inverse") %>% 
    yac_str() %>% 
    as_r()
##      [,1] [,2]
## [1,] "1"  "a" 
## [2,] "0"  "1"

yacas variables

yacas also has variables. They are assigned by :=.

Consider this example:

yac_str("poly := (x-3)*(x+2)")
## [1] "(x-3)*(x+2)"

If the output is not necessary, it can be suppressed by using yac_silent() instead of yac_str():

yac_silent("poly := (x-3)*(x+2)")

We can now list yacas variables (I is the imaginary unit):

yac_str("Variables()")
## [1] "{t,rformBitwiseOps,TeXForm'FuncPrec,I,j,poly,Arow,nd}"
yac_str("Expand(poly)")
## [1] "x^2-x-6"
"poly" %>% y_fn("Expand") %>% yac_str()
## [1] "x^2-x-6"

Sums

Yacas can sum an expression. The syntax is Sum(var, from, to, body) as described in the yacas documentation of Sum(). For example we can sum \(a^k\) for \(k = 0\) to \(k = n\) as follows:

yac_str("Sum(k, 0, n, a^k)")
## [1] "(1-a^(n+1))/(1-a)"

Limits

Yacas can also take the limit of an expression. The syntax is Limit(var, val) expr as described in the yacas documentation of Limit(). For example we can take the limit of \((1+(1/n))^n\) for \(n \to \infty\) as follows:

cmd <- "Limit(n, Infinity) (1+(1/n))^n"
yac_str(cmd)
## [1] "Exp(1)"
yac_expr(cmd)
## expression(exp(1))

This can also be used to illustrate taking derivatives, e.g. the derivative of sine:

yac_str("Limit(h, 0) (Sin(x+h)-Sin(x))/h") 
## [1] "Cos(x)"

Solving equations

Say we want to find roots of poly. First note that the equality in yacas is ==. Then:

cmd <- "Solve(poly == 0, x)"
cmd %>% yac_str()
## [1] "{x==(-2),x==3}"

Note that the default in yacas’s Solve() is to find roots if no equality sign is provided:

cmd <- "Solve(poly, x)"
cmd %>% yac_str()
## [1] "{x==(-2),x==3}"

If we want the solution without the variable name x and the equality symbol ==, we can use Ryacas helper function y_rmvars():

cmd
## [1] "Solve(poly, x)"
cmd %>% y_rmvars() %>% yac_str()
## [1] "{-2,3}"
cmd %>% y_rmvars() %>% yac_expr()
## expression(c(-2, 3))

We can also use y_fn():

"poly == 0" %>% y_fn("Solve", "x")
## [1] "Solve(poly == 0, x)"
"poly" %>% y_fn("Solve", "x") # default is == 0
## [1] "Solve(poly, x)"
"poly == 0" %>% y_fn("Solve", "x") %>% y_rmvars() %>% yac_str()
## [1] "{-2,3}"

Case: gradient and Hessian for a function

Say we have a function:

f <- function(x, y) 2*x^2 + 3*y + y*x^2
f_body <- body(f)
f_body
## 2 * x^2 + 3 * y + y * x^2
f_body_chr <- as.character(as.expression(body(f)))
f_body_chr
## [1] "2 * x^2 + 3 * y + y * x^2"
# or:
# f_body_chr <- "2 * x^2 + 3 * y + y * x^2"
# f <- function(x, y) NULL
# body(f) <- parse(text = f_body_chr, keep.source = FALSE)

The gradient can be found using yacas (D(x) expr is the derivative of expr with respect to x):

cmd_g <- paste0("{ D(x) ", f_body_chr, ", D(y) ", f_body_chr, " }")
cmd_g
## [1] "{ D(x) 2 * x^2 + 3 * y + y * x^2, D(y) 2 * x^2 + 3 * y + y * x^2 }"
g_body <- yac_expr(cmd_g)
g_body
## expression(c(4 * x + y * 2 * x, x^2 + 3))
g <- function(x, y) NULL
body(g) <- g_body
g
## function (x, y) 
## c(4 * x + y * 2 * x, x^2 + 3)
g(2, 4)
## [1] 24  7

The Hessian matrix can also be found using yacas:

cmd_H <- paste0("HessianMatrix(", f_body_chr, ", {x, y})")
cmd_H
## [1] "HessianMatrix(2 * x^2 + 3 * y + y * x^2, {x, y})"
H_body <- yac_expr(cmd_H)
H_body
## expression(rbind(c(2 * y + 4, 2 * x), c(2 * x, 0)))
H <- function(x, y) NULL
body(H) <- H_body
H
## function (x, y) 
## rbind(c(2 * y + 4, 2 * x), c(2 * x, 0))
H(2, 4)
##      [,1] [,2]
## [1,]   12    4
## [2,]    4    0

Ryacas low-level reference

Principle:

Reference: