Ryacas
functionalitylibrary(Ryacas0)
Ryacas
makes the yacas
computer algebra
system available from within R
. The name yacas
is short for “Yet Another Computer Algebra System”. The
yacas
program is developed by Ayal Pinkhuis (who is also
the maintainer) and others, and is available at http://www.yacas.org/ for
various platforms. There is a comprehensive documentation (300+ pages)
of yacas
(also available at http://www.yacas.org/) and the documentation contains
many examples.
R
expressions, yacas
expressions and
Sym
objectsThe Ryacas
package works by sending
`commands'' to
yacaswhich makes the calculations and returns the result to
R`.
There are various different formats of the return value as well
R
expressionsA call to yacas
may be in the form of an R
expression which involves valid R calls, symbols or constants (though
not all valid R
expressions are valid). For example:
<- yacas(expression(Factor(x^2-1))) exp1
The result exp1
is not an expression in the
R
sense but an object of class "yacas"
. To
evaluate the resulting expression numerically, we can do
Eval(exp1, list(x=4))
## [1] 15
yacas
expressionsSome commands are not proper R
expressions. For example,
typing
yacas(expression(D(x)Sin(x)))
produces an error. For such cases we can make a specification using
the yacas
syntax:
yacas("D(x)Sin(x)")
## yacas_expression(cos(x))
Sym
objectsProbably the most elegant way of working with yacas
is
by using Sym
objects. A Sym
object is a
yacas
character string that has the “Sym” class. One can
combine Sym
objects with other Sym
objects as
well as to other R
objects using +
,
-
and other similar R
operators.
The function Sym(x)
coerces an object x
to
a Sym
object by first coercing it to character and then
changing its class to “Sym”:
<- Sym("x") x
Operations on Sym
objects lead to new Sym
objects:
+4 x
## yacas_expression(x + 4)
One can apply sin
, cos
, tan
,
deriv
, Integrate
and other provided functions
to Sym
objects. For example:
Integrate(sin(x), x)
## yacas_expression(-cos(x))
In this way the communication with yacas
is
``tacit’’.
It is important to note the difference between the R
name x
and the symbol "x"
as illustrated
below:
<- Sym("xs")
x x
## yacas_expression(xs)
+4 x
## yacas_expression(xs + 4)
Eval(x+4, list(xs=5))
## [1] 9
The convention in the following is 1) that Sym
objects
match with their names that they end with an ‘s’, e.g.
<- Sym('xs') xs
Algebraic calculations:
yacas(expression((10 + 2) * 5 + 7^7))
## yacas_expression(823603)
yacas(expression(1/14+5/21* (30- 1+ 1/2)))
## yacas_expression(149/21)
#(Sym(10) + Sym(2)) * Sym(5) + Sym(7) ^ Sym(7)
Sym("10 * 2") * 5 + Sym(7) ^ 7
## yacas_expression(823643)
#(Sym(10) + 2) * 5 + Sym(7) ^ 7
#Sym("(10+2)*5 + 7^7")
Sym("1/14 + 5/21 * (30 - 1+1/2)")
## yacas_expression(149/21)
Numerical evaluations:
yacas(expression(N(-12/2)))
## yacas_expression(-6)
Sym("-12/2")
## yacas_expression(-6)
#Eval(Sym("-12/2"))
Symbolic expressions:
yacas(expression(Factor(x^2-1)))
## yacas_expression((x + 1) * (x - 1))
<- expression(x^2 + 2 * x^2)
exp1 <- expression(2 * exp0)
exp2 <- expression(6 * pi * x)
exp3 <- expression((exp1 * (1 - sin(exp3))) / exp2)
exp4 yacas(exp4)
## yacas_expression(3 * (x^2 * (1 - sin(6 * (x * pi))))/(2 * exp0))
Factor(xs^2-1)
## yacas_expression((xs + 1) * (xs - 1))
<- xs^2 + 2 * xs^2
exp1 <- Sym("exp0")
exp0 <- 2 * Sym(exp0)
exp2 <- 6 * Pi * xs
exp3 <- exp1 * (1 - sin(exp3)) / exp2
exp4 exp4
## yacas_expression(3 * (xs^2 * (1 - sin(6 * (xs * pi))))/(2 * exp0))
Combining symbolic and numerical expressions:
yacas(expression(N(Sin(1)^2 + Cos(x)^2)))
## yacas_expression(cos(x)^2 + 0.7080734182)
N(sin(1)^2 + cos(xs)^2)
## yacas_expression(cos(xs)^2 + 0.708073418273571)
Differentiation:
yacas("D(x)Sin(x)")
## yacas_expression(cos(x))
deriv(sin(xs), xs)
## yacas_expression(cos(xs))
Integration:
yacas("Clear(a,b,A,B)")
## yacas_expression(TRUE)
yacas("Integrate(x,a,b)Sin(x)")
## yacas_expression(cos(a) - cos(b))
<- Sym("as")
as <- Sym("bs")
bs Integrate(sin(xs), xs, as, bs)
## yacas_expression(cos(as) - cos(bs))
Expanding polynomials:
yacas("Expand((1+x)^3)")
## yacas_expression(x^3 + 3 * x^2 + 3 * x + 1)
Expand((1+xs)^3)
## yacas_expression(xs^3 + 3 * xs^2 + 3 * xs + 1)
Taylor expansion:
yacas("texp := Taylor(x,0,3) Exp(x)")
## yacas_expression(x + x^2/2 + x^3/6 + 1)
<- Taylor(exp(xs), xs, 0, 3)
texp # Or: texp <- Sym("texp")
Printing the result in nice forms:
yacas("PrettyForm(texp)")
##
## 2 3
## x x
## x + -- + -- + 1
## 2 6
yacas("TeXForm(texp)")
## [1] "x + \\frac{x ^{2}}{2} + \\frac{x ^{3}}{6} + 1"
PrettyForm(texp)
##
## 2 3
## xs xs
## xs + --- + --- + 1
## 2 6
TeXForm(texp)
## [1] "xs + \\frac{xs ^{2}}{2} + \\frac{xs ^{3}}{6} + 1"
yacas
calculationsThe function Set()
and the operator :=
can
both be used to assign values to global variables.
yacas("n := (10 + 2) * 5")
## yacas_expression(60)
yacas("n := n+n")
## yacas_expression(120)
#yacas("Set(z, Cos(a))")
#yacas("z+z")
The same can be achieved with Sym
objects: Consider:
Set(ns, (10 + 2) * 5)
## yacas_expression(60)
Now ns
exists as a variable in yacas
(and
we can make computations on this variable as above). However we have no
handle on this variable in R
. Such a handle is obtained
with
<- Sym("ns") ns
Now the R
variable ns
refers to the
yacas
variable ns
and we can make calculations
directly from R
, e.g:
Set(ns,123)
## yacas_expression(123)
ns
## yacas_expression(123)
^2 ns
## yacas_expression(15129)
Likewise:
<- Sym("as")
as <- Sym("zs")
zs Set(zs, cos(as))
## yacas_expression(cos(as))
+ zs zs
## yacas_expression(2 * cos(as))
Clear a variable binding execute Clear()
:
yacas(expression(n))
## yacas_expression(120)
yacas("Clear(n)")
## yacas_expression(TRUE)
yacas(expression(n))
## yacas_expression(n)
Set(ns, 1)
## yacas_expression(1)
<- Sym("ns")
ns ns
## yacas_expression(1)
Clear(ns)
## yacas_expression(TRUE)
ns
## yacas_expression(ns)
Evaluations are generally exact:
yacas("Exp(0)")
## yacas_expression(1)
yacas("Exp(1)")
## yacas_expression(exp(1))
yacas("Sin(Pi/4)")
## yacas_expression(root(1/2, 2))
yacas("355/113")
## yacas_expression(355/113)
exp(Sym(0))
## yacas_expression(1)
exp(Sym(1))
## yacas_expression(exp(1))
sin(Pi/4)
## yacas_expression(root(1/2, 2))
Sym("355/113")
## yacas_expression(355/113)
To obtain a numerical evaluation (approximation), the
N()
function can be used:
yacas("N(Exp(1))")
## yacas_expression(2.7182818284)
yacas("N(Sin(Pi/4))")
## yacas_expression(0.7071067811)
yacas("N(355/113)")
## yacas_expression(3.1415929203)
N(exp(1))
## yacas_expression(2.71828182845905)
N(sin(Pi/4))
## yacas_expression(0.7071067811)
N(355/113)
## yacas_expression(3.14159292035398)
The N()
function has an optional second argument, the
required precision:
yacas("N(355/133,20)")
## yacas_expression(2.66917293233083)
N("355/113",20)
## yacas_expression(3.14159292035398)
The command Precision(n)
can be used to specify that all
floating point numbers should have a fixed precision of n digits:
yacas("Precision(5)")
## yacas_expression(Precision(5))
yacas("N(355/113)")
## yacas_expression(3.1415929203)
Precision(5)
## yacas_expression(Precision(5))
N("355/113")
## yacas_expression(3.1415929203)
Rational numbers will stay rational as long as the numerator and denominator are integers:
yacas(expression(55/10))
## yacas_expression(11/2)
Sym("55 / 10")
## yacas_expression(11/2)
Some exact manipulations :
yacas("1/14+5/21*(30-(1+1/2)*5^2)")
## yacas_expression(-12/7)
yacas("0+x")
## yacas_expression(x)
yacas("x+1*y")
## yacas_expression(x + y)
yacas("Sin(ArcSin(alpha))+Tan(ArcTan(beta))")
## yacas_expression(alpha + beta)
Sym("1/14+5/21*(1*30-(1+1/2)*5^2)")
## yacas_expression(-12/7)
<- Sym("xs")
xs <- Sym("ys")
ys 0+xs
## yacas_expression(xs)
+1*ys xs
## yacas_expression(xs + ys)
sin(asin(xs))+tan(atan(ys))
## yacas_expression(xs + ys)
The imaginary unit \(i\) is denoted I and complex numbers can be entered as either expressions involving I or explicitly Complex(a,b) for a+ib.
yacas("I^2")
## yacas_expression(-1)
yacas("7+3*I")
## yacas_expression(complex_cartesian(7, 3))
yacas("Conjugate(%)")
## yacas_expression(complex_cartesian(7, -3))
yacas("Exp(3*I)")
## yacas_expression(complex_cartesian(cos(3), sin(3)))
^2 I
## yacas_expression(-1)
7+3*I
## yacas_expression(complex_cartesian(7, 3))
Conjugate(7+3*I)
## yacas_expression(complex_cartesian(7, -3))
exp(3*I)
## yacas_expression(complex_cartesian(cos(3), sin(3)))
%
operatorThe operator %
automatically recalls the result from the
previous line.
yacas("(1+x)^3")
## yacas_expression((x + 1)^3)
yacas("%")
## yacas_expression((x + 1)^3)
yacas("z:= %")
## yacas_expression((x + 1)^3)
1+x)^3 (
## yacas_expression((xs + 1)^3)
<- Sym("%")
zs zs
## yacas_expression((xs + 1)^3)
PrettyForm
, PrettyPrint
and
TeXForm
There are different ways of displaying the output.
The (standard) yacas form is:
yacas("A:={{a,b},{c,d}}")
## Yacas matrix:
## [,1] [,2]
## [1,] a b
## [2,] c d
yacas("B:= (1+x)^2+k^3")
## yacas_expression((x + 1)^2 + k^3)
yacas("A")
## Yacas matrix:
## [,1] [,2]
## [1,] a b
## [2,] c d
yacas("B")
## yacas_expression((x + 1)^2 + k^3)
<- Sym("as"); bs <- Sym("bs"); cs <- Sym("cs"); ds <- Sym("ds")
as <- List(List(as,bs), List(cs,ds))
A <- Sym("ks")
ks <- (1+xs)^2+ks^3
B A
## Yacas matrix:
## [,1] [,2]
## [1,] as bs
## [2,] cs ds
B
## yacas_expression((xs + 1)^2 + ks^3)
The Pretty form is:
yacas("PrettyForm(A)")
##
## / \
## | ( a ) ( b ) |
## | |
## | ( c ) ( d ) |
## \ /
yacas("PrettyForm(B)")
##
## 2 3
## ( x + 1 ) + k
PrettyForm(A)
##
## / \
## | ( as ) ( bs ) |
## | |
## | ( cs ) ( ds ) |
## \ /
PrettyForm(B)
##
## 2 3
## ( xs + 1 ) + ks
The output can be displayed in TeX form:
yacas("TeXForm(B)")
## [1] "\\left( x + 1\\right) ^{2} + k ^{3}"
TeXForm(B)
## [1] "\\left( xs + 1\\right) ^{2} + ks ^{3}"
yacas("40!")
## yacas_expression(8.15915283247898e+47)
Factorial(40)
## yacas_expression(Factorial(40))
Expand \(\exp(x)\) in three terms around \(0\) and \(a\):
yacas("Taylor(x,0,3) Exp(x)")
## yacas_expression(x + x^2/2 + x^3/6 + 1)
yacas("Taylor(x,a,3) Exp(x)")
## yacas_expression(exp(a) + exp(a) * (x - a) + (x - a)^2 * exp(a)/2 + (x - a)^3 * exp(a)/6)
<- Sym("xs")
xs Taylor(exp(xs),xs,0,3)
## yacas_expression(xs + xs^2/2 + xs^3/6 + 1)
<- Sym("as")
as Taylor(exp(x),x,as,3)
## yacas_expression(exp(as) + exp(as) * (xs - as) + (xs - as)^2 * exp(as)/2 + (xs - as)^3 * exp(as)/6)
The InverseTaylor()
function builds the Taylor series
expansion of the inverse of an expression. For example, the Taylor
expansion in two terms of the inverse of \(\exp(x)\) around \(x=0\) (which is the Taylor expansion of
\(Ln(y)\) around \(y=1\)):
yacas("InverseTaylor(x,0,2)Exp(x)")
## yacas_expression(x - 1 - (x - 1)^2/2)
yacas("Taylor(y,1,2)Ln(y)")
## yacas_expression(y - 1 - (y - 1)^2/2)
InverseTaylor(exp(xs),xs,0,2)
## yacas_expression(xs + xs^2/2 + 1)
Taylor(log(ys),ys,1,2)
## yacas_expression(ys - 1 - (ys - 1)^2/2)
Solve equations symbolically with the Solve()
function:
yacas("Solve(x/(1+x) == a, x)")
## Yacas vector:
## [1] x == a/(1 - a)
yacas("Solve(x^2+x == 0, x)")
## Yacas vector:
## [1] x == 0 x == -1
Solve(xs/(1+xs) == as, xs)
## Yacas vector:
## [1] xs == as/(1 - as)
Solve(xs^2+xs == 0, xs)
## Yacas vector:
## [1] xs == 0 xs == -1
Solve(List(xs^2+ys^2==6, xs-ys==3), List(xs,ys))
## Yacas matrix:
## [,1] [,2]
## [1,] xs == root(6 - ys^2, 2) ys == ys
# FIXME
#Solve(List(mean==(xs/(xs+ys)), variance==((xs*ys)/(((xs+ys)^2) * (xs+ys+1)))),
# List(xs,ys))
(Note the use of the == operator, which does not evaluate to anything, to denote an “equation” object.)
To solve an equation (in one variable) like \(sin(x)-exp(x)=0\) numerically taking \(0.5\) as initial guess and an accuracy of \(0.0001\) do:
yacas("Newton(Sin(x)-Exp(x),x, 0.5, 0.0001)")
## yacas_expression(-3.1830630118)
Newton(sin(xs)-exp(xs),xs, 0.5, 0.0001)
## yacas_expression(-3.1830630118)
yacas(expression(Expand((1+x)^3)))
## yacas_expression(x^3 + 3 * x^2 + 3 * x + 1)
Expand((x+1)^3)
## yacas_expression(xs^3 + 3 * xs^2 + 3 * xs + 1)
The function Simplify() attempts to reduce an expression to a simpler form.
<- Sym("y")
y yacas("(x+y)^3-(x-y)^3")
## yacas_expression((x + y)^3 - (x - y)^3)
yacas("Simplify(%)")
## yacas_expression(6 * (x^2 * y) + 2 * y^3)
+y)^3-(x-y)^3 (x
## yacas_expression((xs + y)^3 - (xs - y)^3)
Simplify("%")
## yacas_expression(6 * (xs^2 * y) + 2 * y^3)
Analytical derivatives of functions can be evaluated with the
D()
and deriv()
functions:
yacas("D(x) Sin(x)")
## yacas_expression(cos(x))
deriv(sin(xs), xs)
## yacas_expression(cos(xs))
These functions also accepts an argument specifying how often the derivative has to be taken, e.g:
yacas("D(x,4)Sin(x)")
## yacas_expression(sin(x))
deriv(sin(xs),xs,4)
## yacas_expression(sin(xs))
@ <<echo=F,results=hide>>= yacas(“Clear(a,b,A,B)”)
```r
#yacas("Integrate(x,a,b)Sin(x)")
#yacas("Integrate(x,a,b)Ln(x)+x")
#yacas("Integrate(x)1/(x^2-1)")
yacas("Integrate(x)Sin(a*x)^2*Cos(b*x)")
## yacas_expression((2 * sin(b * x)/b - (sin(-2 * (x * a) - b * x)/(-2 * a - b) + sin(-2 * (x * a) + b * x)/(-2 * a + b)))/4)
#Integrate(sin(x),x,a,b)
#Integrate(log(x),x,a,b)
#Integrate(1/(x^2-1),x)
<- Sym("a")
a <- Sym("b")
b Integrate(sin(a*x)^2*cos(b*x),x)
## yacas_expression((2 * sin(b * xs)/b - (sin(-2 * (xs * a) - b * xs)/(-2 * a - b) + sin(-2 * (xs * a) + b * xs)/(-2 * a + b)))/4)
#yacas("Limit(x,0)Sin(x)/x")
yacas("Limit(n,Infinity)(1+(1/n))^n")
## yacas_expression(exp(1))
yacas("Limit(h,0) (Sin(x+h)-Sin(x))/h")
## yacas_expression(cos(x))
#Limit(sin(x)/x,x,0)
<- Sym("ns")
ns Limit((1+(1/ns))^ns,ns,Infinity)
## yacas_expression(exp(1))
<- Sym("hs")
hs Limit((sin(xs+hs)-sin(xs))/hs,hs,0)
## yacas_expression(cos(xs))
yacas("Subst(x,Cos(a))x+x")
## yacas_expression(2 * cos(a))
Subst(xs+xs,xs,cos(as))
## yacas_expression(2 * cos(as))
yacas("OdeSolve(y''==4*y)")
## yacas_expression(C545 * exp(2 * x) + C549 * exp(-2 * x))
yacas("OdeSolve(y'==8*y)")
## yacas_expression(C579 * exp(8 * x))
yacas("E4:={ {u1,u1,0},{u1,0,u2},{0,u2,0} }")
## Yacas matrix:
## [,1] [,2] [,3]
## [1,] u1 u1 0
## [2,] u1 0 u2
## [3,] 0 u2 0
yacas("PrettyForm(E4)")
##
## / \
## | ( u1 ) ( u1 ) ( 0 ) |
## | |
## | ( u1 ) ( 0 ) ( u2 ) |
## | |
## | ( 0 ) ( u2 ) ( 0 ) |
## \ /
<- Sym("u1")
u1 <- Sym("u2")
u2 <- List(List(u1, u1, 0), List(u1, 0, u2), List(0, u2, 0))
E4 PrettyForm(E4)
##
## / \
## | ( u1 ) ( u1 ) ( 0 ) |
## | |
## | ( u1 ) ( 0 ) ( u2 ) |
## | |
## | ( 0 ) ( u2 ) ( 0 ) |
## \ /
yacas("E4i := Inverse(E4)")
## Yacas matrix:
## [,1] [,2]
## [1,] u1 * u2^2/(u1^2 * u2^2) 1/u1 - u1 * u2^2/(u2^2 * u1^2)
## [2,] 1/u1 - u1 * u2^2/(u2^2 * u1^2) u1 * u2^2/(u2^2 * u1^2) - 1/u1
## [3,] -(u1 * u2/u1)/u2^2 u1 * u2/(u1 * u2^2)
## [,3]
## [1,] -(u1 * u2/u1)/u2^2
## [2,] u2 * u1/(u2^2 * u1)
## [3,] u1/u2^2
yacas("Simplify(E4i)")
## Yacas matrix:
## [,1] [,2] [,3]
## [1,] 1/u1 0 -1/u2
## [2,] 0 0 1/u2
## [3,] -1/u2 1/u2 u1/u2^2
yacas("PrettyForm(Simplify(E4i))")
##
## / \
## | / 1 \ ( 0 ) / -1 \ |
## | | -- | | -- | |
## | \ u1 / \ u2 / |
## | |
## | ( 0 ) ( 0 ) / 1 \ |
## | | -- | |
## | \ u2 / |
## | |
## | / -1 \ / 1 \ / u1 \ |
## | | -- | | -- | | --- | |
## | \ u2 / \ u2 / | 2 | |
## | \ u2 / |
## \ /
<- Inverse(E4)
E4i Simplify(E4i)
## Yacas matrix:
## [,1] [,2] [,3]
## [1,] 1/u1 0 -1/u2
## [2,] 0 0 1/u2
## [3,] -1/u2 1/u2 u1/u2^2
PrettyForm(Simplify(E4i))
##
## / \
## | / 1 \ ( 0 ) / -1 \ |
## | | -- | | -- | |
## | \ u1 / \ u2 / |
## | |
## | ( 0 ) ( 0 ) / 1 \ |
## | | -- | |
## | \ u2 / |
## | |
## | / -1 \ / 1 \ / u1 \ |
## | | -- | | -- | | --- | |
## | \ u2 / \ u2 / | 2 | |
## | \ u2 / |
## \ /
yacas("Determinant(E4)")
## yacas_expression(-(u1 * u2^2))
yacas("Determinant(E4i)")
## yacas_expression((u1 * u2^2/(u2^2 * u1^2) - 1/u1) * (u1 * u2^2) * u1/(u1^2 * u2^4) - u1 * u2^2 * (u2 * u1) * (u1 * u2)/(u1 * u2^2 * (u2^2 * u1 * (u1^2 * u2^2))) - (1/u1 - u1 * u2^2/(u2^2 * u1^2)) * (u1 * u2)^2/(u2^4 * u1^2) - (1/u1 - u1 * u2^2/(u2^2 * u1^2))^2 * u1/u2^2 - (1/u1 - u1 * u2^2/(u2^2 * u1^2)) * (u2 * u1) * (u1 * u2)/(u2^4 * u1^2) - (u1 * u2^2/(u2^2 * u1^2) - 1/u1) * (u1 * u2)^2/(u1^2 * u2^4))
yacas("Simplify(E4i)")
## Yacas matrix:
## [,1] [,2] [,3]
## [1,] 1/u1 0 -1/u2
## [2,] 0 0 1/u2
## [3,] -1/u2 1/u2 u1/u2^2
yacas("Simplify(Determinant(E4i))")
## yacas_expression(-1/(u1 * u2^2))
determinant(E4)
## yacas_expression(-(u1 * u2^2))
determinant(E4i)
## yacas_expression((u1 * u2^2/(u2^2 * u1^2) - 1/u1) * (u1 * u2^2) * u1/(u1^2 * u2^4) - u1 * u2^2 * (u2 * u1) * (u1 * u2)/(u1 * u2^2 * (u2^2 * u1 * (u1^2 * u2^2))) - (1/u1 - u1 * u2^2/(u2^2 * u1^2)) * (u1 * u2)^2/(u2^4 * u1^2) - (1/u1 - u1 * u2^2/(u2^2 * u1^2))^2 * u1/u2^2 - (1/u1 - u1 * u2^2/(u2^2 * u1^2)) * (u2 * u1) * (u1 * u2)/(u2^4 * u1^2) - (u1 * u2^2/(u2^2 * u1^2) - 1/u1) * (u1 * u2)^2/(u1^2 * u2^4))
Simplify(E4i)
## Yacas matrix:
## [,1] [,2] [,3]
## [1,] 1/u1 0 -1/u2
## [2,] 0 0 1/u2
## [3,] -1/u2 1/u2 u1/u2^2
Simplify(determinant(E4i))
## yacas_expression(-1/(u1 * u2^2))
Note that the value returned by yacas
can be of
different types:
yacas(expression(Factor(x^2-1)),retclass='unquote')
## *" ("+" (x ,1 ),"- (x ,1 ))
yacas(expression(Factor(x^2-1)),retclass='character')
## "*" ("+" (x ,1 ),"-" (x ,1 ))
Set output width:
get_output_width()
## [1] 60
set_output_width(120)
get_output_width()
## [1] 120