4. Defining Functions. Substitutions.  Delayed 	  	Assignments. Transforming Expressions.

2016-04-15

$Version

4.1 User Defined Functions

4.1.1	Definition of a function

Variables with an underscore ("_")  at the left hand sinde are dummy variables. Other variables are global.

Clear[f, x,y,c];
f[x_,y_] = x^2 + y^ 2 - c

In the example above  x  and  y  are dummy variables,  c  is global. It is important that no assignements have been made for the dummy variables used in the definition of the function. This is explained in more detail below.

f[a,b]

f[2,5]

c = 13.3

f[2,5]

f[a,b]

Clear[c]
f[a,b]

f[x_] = x^3 + x - 1

sf = NSolve[f[x] == 0,x]

Now the consequences are shown, which result from a variable with a previous assignement (here  x) used in the definition of a function.

x = 37

f[x_] = x^3 + x - 1

f[d]

4.1.2	Properties of Functions: Selfcalls

A function may call itself  in Mathematica. This is a feature forbidden in some other programming
linguages. This may be used to calculate an expression or a function defined by a recurrence which calls this same function. For example, the Legendre polynomials  Pn(x)  fullfilling the recurrence 
                           (n + 1) Pn+1(x) - (2 n + 1) x Pn(x)  +  n Pn-1(x)  =  0

may be generated in the following way:

Clear[p,x]
p[0,x_] = 1;  p[1,x_] = x;
p[(n_?NonNegative)?IntegerQ,x_]:=
p[n,x]=((2n-1)/n x) p[n-1,x] - p[n-2,x] (n-1)/n;

p[5,z]//Apart

LegendreP[5,z]//Apart

p[.5,z]

LegendreP[.5,z]

LegendreP[-6,z]

Though the recurrence is valid for any  (real or complex)  value   n , the above program gives correct results for a postive integer  n, only.  This is safeguarded by the checks attached to the argument   n   at the left hand side.

4.1.3	Properties of Functions: Dynamic Programming

Recurrences can be defined in two ways. The second one of  those given below is much faster and is called dynamic programming. The recurrence below gives the Fibonacci numbers; in the sequence each number is the sum of the two preceeding numbers. The first two numbers are 1.

Clear[fib];  fib[0]=fib[1]=1;
fib[n_]:=fib[n-1]+fib[n-2];
nt = Table[{k,fib[k]}//Timing,{k,2,25}]; 
nt[[Range[1,14]]]

The first number in each sublist is the computing time needed to compute the Fibonacci
number given as the very last number in each sublist.

tp=Table[{nt[[k,2,1]],nt[[k,1]]},{k,Length[nt]}];
ListPlot[tp,PlotRange->All,AxesLabel->{"n","time [s]"},PlotLabel->"time for computing the n-th \nFibonacci number\n\n"]

Clear[fib]; fib[0]=fib[1]=1;
fib[n_]:=fib[n]=fib[n-1]+fib[n-2] ;
dt = Table[{k,fib[k]},{k,2,200}]//Timing ;  dt[[1]]

The preceeding block of expressions presents the same recurrence in dynamic programming.In dynamic programming it takes almost the same time to compute the first 200 numbers as is needed to compute the 14-th Fibonacci number in the first approach ! The latter times increase steeply with n; e.g. for n = 41, it amounts to 13324 s = 3.7 h  !   Mathematica has also a programme for computing Fibonacci numbers:

dm=Table[{k,Fibonacci[k]},{k,2,200}]//Timing; dm[[1]]

4.1.4 Functions with two list of arguments

Clear[a,b,c,x]
expr=a*x+b*x^2
D[expr,x]

%/.{a->1,b->2,x->c}

Here is another method to get the same result;  it uses a function definition for the original expression. By using separate square brackets for the parameters, (a,b),  and the variable,x, one can use the standard Mathematica 
derivative notation.

Clear[f,a,b,x]
f[a_,b_][x_]:=a*x+b*x^2
f[1,2]'[c]

4.1.5	Substitutions in  Function Definitions

In the function below any addition is replaced by a multiplication. To understand this one should know a litte 
how Mathematica  expressions are stored (Chap.20).

Clear[g]
g[x_] := x /. Plus -> Times

g[a + b + c]

4.1.6	Conditions in  Function Definitions

Conditions may be introduced into function definitions by use of   /; another method may use branching
commands ( If[], Which[],...), s. sect. 7.4.

p[x_]:=(x^2)/;x>0;p[x_]:=-(x^2)/;x≤0

{p[1],p[.5], p[0], p[-.5],p[-1]}

There is still another way: Introducing conditions into the argument of the function:

Clear[f,n]

f[(n_?Positive)?IntegerQ] = n!
	f[n_] := Print["f expects a positive integer argument"]

f[3]

f[-4]

f[2.2]

4.2   Built-In  Functions.

Mathematica contains numerous programs for the evaluation of elementary or special functions. These contain several rules for analytic computations. These yield numeric values for numeric arguments; these may be taken from wide domains.  A table of the names for the most important functions is given in the file MathFunctionsinMMA.pdf on the website.

4.3  Pure Functions

The definition of a functions in the first paragraph involved dummy variable denoted by underscores  in the list of arguments. These variables are not called by their names; the names only serve as symbols to denote a position or variable. Another way to achieve such an aim is the pure function. Here we give a short introduction; the subject will be treated in more detail insect. 21.3; in particular the use of  pure functions with the operators  Map[] and Apply[].

body &      a pure function in which arguments are specified as  #  or  #1,#2,#3,... etc.

The ampersand     &  is obligatory in this way  of defining a pure function. Examples are:

x //#^2 &

a //#^2 &

Pi//N[#,22]&

The pure functions are applied in postfix form (sect.20.2) to expressions preceeding them. Another important application of pure functions is in searching for or selecting elements in lists,  see sect. 5.8

4.4 Substitutions  (Rule)

Substitutions must be prescribed according to the following rules:

exp   /. 	var -> value

	 exp   /.	list

			list = {var1 -> val1, var2 -> val2, ... }

A variable at the left of an arrow occuring in  exp  is assigned the value given at the right of the arrow; so a new expression is obtained. The old expression remains unchanged.

h = z

h1 = h /. z -> 44

h

h1

Clear[x,y,f]
f[x_,y_] = x^2 + y^ 2 - c

h = f[x,y]

h1 = h /. c -> 13.3

h

h1

h2 = h1 /. {x -> a, y -> b}

x = 4
y = 5

h

h1

h2

f[p,q]

f[x,y]

4.4.1	Iterated Substitutions

Two slashes in the substitution command incite Mathematica to apply the substitution several times:

Clear[a,b,c,d,e,f,u,v,w, x,y]

u = (a + b)^2; v = (c + d)^2;  w = (e + f)^2;

uv = Expand[u + v + w]

uv /. x_^2 + 2 x_ y_ + y_^2 -> (x + y)^2

uv //. x_^2 + 2 x_ y_ + y_^2 -> (x + y)^2

4.4.2	 Substitutions Restricted by Conditions

Substitutions may be restricted by the command  /; as it was   done  in the definitinion of functions.

5!

5!!   (* 5 x 3 x 1 *)

Clear[m,n,k]

Γ(m+ 3) /. Γ(n_+ k_ /; OddQ[2k]) -> (((2 k+2 n-2)!! (Sqrt[π]))/(2^(k-(1/2)+n)))

Γ(m+(5/2)) /. Γ(n_+ k_ /; OddQ[2k]) -> (((2 k+2 n-2)!! (Sqrt[π]))/(2^(k-(1/2)+n)))

% /. m -> 4

%% /. m -> 4.3

4.4.3	 Delayed Substitutions (Rule Delayed)

lhs :> rhs  	is a substitution that transforms lhs to rhs, evaluating rhs only after the 			rule is applied.
:>   is generated by the colon (:)  followed by the greater-than sign (>)

su = x :> t

f = x^2

f /. su

t = 5;
f /. su

Clear[t]
f /. su

4.5 Immediate and Delayed Definitions (Assignments)

lhs = rhs  	Immediate assignment: rhs is evaluated when the assignment is made.
		rhs is intended to be the final value of the name lhs.
		
lhs := rhs 	Delayed assignment: rhs is evaluated each time the value of  lhs is requested.
		 rhs gives a "command" or "program" to be executed whenever one asks for 
		 the value of   lhs.

Clear[p, s, x, y]

x = 4

s = x^2

p := x^2

x = 5

Print[s]

Print[p]

Clear[x]; ex[x_] := Expand[(1 + x)^2]

?ex

iex[x_] = Expand[(1 + x)^2]

?iex

rd = ex[y + 2]

ri = iex[y + 2]

rd - ri

4.5.1	 Further applications of assignments

The transcendental equation

f1(x,a) = x - exp(-a x) = 0

has a single real root  (Subscript[x, 0])  =  (Subscript[x, 0])(a), which is a function of the parameter  a.  The program below calculates this root for a given value  a  by  a numeric method.

f1[x_,a_] := - Exp[- a x] + x ; 
g0[a_]    := If[a > 0, 1/a, 1/(1 + a) ];
g1[a_]    := FindRoot[f1[x,a], {x,g0[a]} ];

Plot[x/.g1[a],{a,-0.35`,15}]

f1[x_,a_] = - Exp[- a x] + x 
g0[a_]    = If[a > 0, 1/a, 1/(1 + a) ]
g1[a_]    = FindRoot[f1[x,a], {x,g0[a]} ]

Plot[x /. g1[a], {a,-.35,15}]

The result is same figure as above. So one can work with immediate or delayed assignments. Using the latter avoids the error messages and other unnecessary output.

4.6    Transforming Alebraic Expressions

Expand[expr]		multiply out products and powers, writing the result as a sum of  terms  
				            
ExpandAll[expr] 	apply  Expand[]  everywhere				   

Factor[expr]  		reduce to a product of factors	   	   
	
Together[expr]   	put all terms over a common denominator 

Apart[expr]    	separate into terms with simple denominators	
Apart[expr,var]    	partial fraction decomposition w.r.t. variable   var	   

Cancel[expr]     	cancel common factors between  numerators and denominators         
                  											   	   
Simplify[expr]		try  a sequence of algebraic transformations and give the smallest form 
				of    expr   found	
			
Simplify[expr,assump]	does simplification using assumptions assump 
	                   		 	
FullSimplify[expr]	may lead to still simpler expressions than  Simplify[],  but is sometimes 
				very time-consuming

Mathematica  contains these commands in the menue "Algebraic Manipulations". 
This is called via the menues "File", "Palettes".

f = x^7 - a^7

g = Factor[f]

h = x - a

k = f/g

Simplify[k]

Cancel[k]

k = f/h

Simplify[k]

Cancel[k]

e = (x - 1)^2 (2 + x) /( (1 + x) (x - 3)^2 )

Expand[e]

ExpandAll[e]

et = Together[%]

ae = Apart[%]

Factor[%]

s = Simplify[%]

n = Numerator[s]

d = Denominator[s]

d = Factor[d]

n/d

e == %

e == ae

e == et

The expression on the rhs was obtained from that on the lhs. But Mathematica does not take note of their identity. The equality must be transformed. One may start with Together[] or ExpandAll[]; if these are not successful one may continue with the more powerful but also more time-consuming operators Simplify[] or even FullSimplify[].

e == et // ExpandAll

e == ae // ExpandAll

e == ae // Simplify

Clear[x,y]
f = E^Abs[x - y]

The derivative of this function exists. But Mathematica cannot handle it directly:

D[f,x]

The above derivative of f is not usable. So one must distinguish the two cases; this is accomplished by the option for assumptions in Simplify[] :

g = Simplify[f, x > y]
D[g,x]

g = Simplify[f, x < y]
D[g,x]

Collect[expr, x] 		group together powers of  x

This command may be used to order complicated expressions according to the various variables. it may be necessary to procedd step by step:

f = 4 x + 6 y + 10 z

g = Collect[f^3, x]

h = Collect[f^3, y]

PowerExpand[expr]     transform (x y)^p  into x^p y^p , etc.

PowerExpand[]  must be used with great care. A sloppy use of it may give wrong results, in particular for complex values or variables.

f = Sqrt[x y]

g = PowerExpand[f]

f - g

Simplify[%]

ExpandAll[%%]

PowerExpand[%%%]

PowerExpand[ Sqrt[- x y]]

Tan[ArcTan[x]]

ArcTan[Tan[x]]

PowerExpand[%]

FullSimplify[%%]

4.6.1  Treating complex expressions

ComplexExpand[expr]   perform expansions assuming that all variables are real  

ComplexExpand[expr, opt] as above but with opt steering the output.   Options: 

TargetFunctions -> {Re, Im}   Rectangular coordinates in the complex plane

TargetFunctions -> {Abs, Ar}  Polar coordinates in the complex plane

In order to save space, input and output are displayed in the same line below:

Input 							Output

z = x + I y                                                                   	       x+i y

Re[z^3]                                                                         Re[(((x+i y))^3)]

Im[z^3]                                                                         Im[(((x+i y))^3)]

rg = ComplexExpand[Re[z^3]]		       x3 - 3 x y2

ig = ComplexExpand[Im[z^3]]			       3 x2 y - y3

rg + I ig                                                                          (x^3)-3 x (y^2)+i (3 (x^2) y-(y^3))

ComplexExpand[z^3]                                               (x^3)-3 x (y^2)+i (3 (x^2) y-(y^3))

z = x + I y

ComplexExpand[1/z^3]

ComplexExpand[Sin[x + I y]]

4.6.2	Further Examples of Transformations and Simplifications

Below is a partical fraction decomposition w.r.t. the variable w, where square root expressions are involved

Clear[n,z1,z2,z3,w]

Apart[((Sqrt[-n+w])/((Sqrt[w]) (Sqrt[w-z1]) (Sqrt[w-z2]) (Sqrt[w-z3]))),w]

In inserting complicated alebraic expressions into an algebraic equation one may encounter
difficulties to simplify these

Clear[u]

f = -2 u+2 (u^3)+2 ε-2 (u^2) ε

df = D[f,u]

so = Solve[df == 0,u]

f0 = f /. so

Expand[f0]

Simplify[f0]

FullSimplify[f0]

FullSimplify[(Sqrt[2 (Sqrt[6])+5])]

Simplify[Sqrt[x^2]]

Simplify[Sqrt[x^2], x > 0]

Simplify[x^2 > 3, x > 2]

Simplify[(m^n)∈Integers,{m,n}∈Integers&&m>0&&n>0]

Simplify[a/b > 0, a > 0 && b > 0]

Simplify[Sqrt[b^2], a*b > 0 && a > 0]

Integrate[ Sin[a x] Cosh[b x]/Sinh[x], {x,0,Infinity}]

(∫_(-∞)^∞)((Sin[a x] Cosh[b x] )/(Sinh[x]))\[DifferentialD]x

Gamma[x] Gamma[1 - x]

FunctionExpand[%]

FunctionExpand[BesselJ[n, I x]]

FunctionExpand[BesselY[n, I x]]

Hypergeometric2F1[1/2,1/2,3/2, Sin[z]^2]

PowerExpand[%]

Hypergeometric2F1[1/2,1,3/2, z^2]

PowerExpand[%]

Clear[n,z,t]

Hypergeometric2F1[-n/2,-(n - 1)/2,1/2, z^2/t^2]

PowerExpand[%]

Hypergeometric2F1[1-n,1,2, -z/t]

PowerExpand[%]

ExpandAll[%]

4.6.3	Treating complex expressions with the operators Simplify[], FunctionExpand[]  may lead to 
	wrong results:

expr = π x/(x + 1 - 2 (-1)^(1/3) + I Sqrt[3])

expr /. x -> 0

Simplify[expr]

Clear[z]
s1 = Hypergeometric2F1[1/2, 1, 2, 4 z (1 - z)]

s2 = FunctionExpand[%]

s3 = PowerExpand[s1]

s4 = Cancel[s1]

{s1, s2, s3, s4} /. z -> 1/4

Hypergeometric2F1[1/2, 1, 2, 4/4 (1 - 1/4)]

4.6.4	Collecting variables with non-integer exponents

expr = Expand[Sum[(-b + a n )x^(n + 0.12/n),{n,3}]]

Collect[expr,x]

Collect[] only works with integer exponents. A way to perform the Collect[] in expressions with non-integer exponents  is to allow a pattern variable for the exponent. This will in effect create a separate "variable" for each distinct power, and this suffices to do what one wants in this particular example.

Collect[expr,(x^(_.))]

(x^(_.))  is generated in the following way: Type x ; use  (■^□) contained in the menue “Palettes” “Other-Basic Typesetting” or “Other-Basic Math Input” ; put the curser into the empty superscript square and type  the keys “underscore” and then “point”.

FullForm[(x^(_.))]

4.6.5	Collecting logarithms

Logarithms can be combinend with the help of the following two commands:

Simplify[Log[x] + Log[y]]

FullSimplify[Log[x] + Log[y]]

Define the following function:

CollectLogs[xx_]:=Log[Simplify[E^xx]]

CollectLogs[Log[x]+Log[y]]

CollectLogs[Log[x]-Log[y]]

or using Simplify[] , provided the two arguments are real and have the same sign:

Simplify[Log[a]-Log[b],Element[{a,b},Reals]&&a>0&&b>0]

Simplify[Log[a]-Log[b],{b>0,a>0}]

Simplify[Log[a]-Log[b],{b<0,a<0}]

4.7    Treating Expressions Containing Trigonometric,
          Hyperpolic Functions and Exponentials

TrigExpand[expr]  	expands out  trigonometric (and/or hyperbolic) functions  in expr

TrigFactor[expr]	  	factors trigonometric functions  in expr

TrigFactorList[expr]		factors trigonometric functions  in expr, yielding a list of lists 								containing trigonometric monomials and exponents.

TrigReduce[expr]	 rewrites products and powers of trigonometric functions in  expr 
			in terms of trigonometric functions with combined  arguments. 

ExpToTrig[expr]	 	converts exponentials in expr to trigonometric functions,
				works also on hyperbolic functions

TrigToExp[expr]	 	converts trigonometric function in expr to exponentials, 
				works also on hyperbolic functions

f = Sin[x]^3  Cos[2 x]

g =TrigExpand[f]

Expand[ g /. Cos[x] -> (1 - Sin[x]^2)^(1/2) ]

TrigFactor[f]

r = TrigReduce[f]

TrigExpand[r]

t = TrigToExp[f]

Expand[t]

ExpToTrig[%]

f

f == %% //Simplify

2+Cos[2 x]+Cos[2 y]+Cos[2 (x+y)]

%/.Cos[2x_]:>2Cos[x]^2-1

4.8	Exercises

4.1          Define the function   f(x,n)=(x^n).    Evaluate f(2,1),f(3,2),f(4,7),f(y,k).

4.2          Transform the expression   f=sin(k1 x) sin(k2 y) into  g by replacing k1  with  a   and  k2  with  b.

4.3	Decompose the following expressions into partial fractions; thereafter put them over common 
            denominator and expand completely. At the end  simplify all these expression as much as possible.

1)         (((x^2)+1)/((x-2) ((((x^2)+1))^2)))  ;                  2 )    (((x^3)+3 (x^2)-4 x+3)/(((x^2)-1) ((((x^2)+1))^2)))     ;       3)      (x/((x^4)-1))

4.4           Plot the function    f(x)  =  sin(x)/(1+x^2)    in the intervall     (0,π )    and determine the 
	maximum  xm    and    f(xm).

4.5          Assuming that x and y are real,compute the real and imaginary parts 
	of the following expressions:

1)      (((x+I y))^5);     2)     cos(x + i y) ;      3)        (((x + i y))^2)sin(x + i y).

4.6	Transform the following expressions into Fourier sums and into pure powers of  sin x   	and  cos x.   These results are not unique in view of the relation   (sin^2)x  +  (cos^2)x  = 1.  	In addition some transformations must be imposed by presenting some trigonometric 	relations as substitutions.

1)   cos(4 x)  (sin^5)x,	2)   (sin^2)(2x) + (cos^2)(2x),	3)  sin(3x) cos(5 x) (cos^2)x;

4.7	A series circuit consists of a resistor with resistance  R, a capacitor with capacitance  
	C  and a coil of inductance L. For a given angular frequency  ω  the impedance of this 	circuit is: Z =  R + i ω L + 1/(i ω C).
	Compute the admittance Y = 1/Z ; decompose it into the real and the imaginary part.

4.8	Spherical Bessel functions  z(m,x) are proportional to Bessel functions with half odd 	integer order:
			(Subscript[j, m])(x) := (Sqrt[(π/(2x))]) (Subscript[J, m + 1/2])(x),  ym(x) := (Sqrt[(π/(2x))]) (Subscript[Y, m + 1/2])(x) ;
	
	they fulfil the following recurrence relations :
	
	z(m+1,x) -  (2m + 1)/x z(m,x) + z(m-1,x)  =  0;	
	j(0,x) = y(-1,x) = sin(x)/x;  j(-1,x) = - y(0,x) = cos(x)/x
	
	Define a function, which computes  j(m,x)  or  y(m,x) for arbitrary natural   m. 
	So z is either j or y. Compute the first few  j's and  y's.

4.9	Define the Heaviside step functionwithout a branching command.

4.10	Get a numeric value of  (e^2)  to 31 decimal places with one postfix command,
	which uses only e as input.

4.11        Expand the sum   (r+s+t+v)    and simplify it again:
	r=(((a+b))^4), s=(((c+d))^4), t =(((e+f))^4), u =(((x+y))^4).

4.12	Find the roots of the following polynomial and verify that they fulfil the corresponding 
	equation p(u) = 0.		p(u) = -2 u+2 (u^3)+3 ε-2 (u^2) ε

4.13	Get the simplest expressions for the following hypergeometric series:
	1)	Hypergeometric2F1[1/2,1/2,3/2, Sin[z]^2]
	2)	Hypergeometric2F1[1,1,3/2, Sin[z]^2]
	3)	Hypergeometric2F1[n/2,-n/2,1/2, Sin[z]^2]
	4)	Hypergeometric2F1[1/2,1,2, 4 z (1 - z)]

4.14	Compute analytical expressions (polynomials in x)  for the first 7 Chebishev 		polynomials from the defining relation   T(x)  =  cos(n arccos(x)), (n = 0,1,2,...,6) .