# Essays/Extremal Arguments

We consider the values of some algebraic and transcendental functions on extremal arguments. The arguments and results are represented as 64-bit IEEE floating point numbers, with the letters x and y denoting real numbers while the letter z denoting complex numbers.

## Utilities

```sin  =: 1&o.
cos  =: 2&o.
sinh =: 5&o.
cosh =: 6&o.

hex  =: (4+#@\$) }. 2&(3!:3)
unhex=: (3!:2) @ ((}: 2 (3!:3) 0.5)&,) " 1
```

The smallest positive IEEE number is 0000000000000001 in hex and the largest is 7fefffffffffffff , respectively  2^_1074 and  (2-2^_52)*2^1023 .

```   unhex (15\$'0'),'1'
4.94066e_324
2^_1074
4.94066e_324
unhex '7fe',13\$'f'
1.79769e308
(2-2^_52)*2^1023
1.79769e308

e=: 2^_1074
E=: (2-2^_52)*2^1023

hex e,E
0000000000000001
7fefffffffffffff
```

## sin and cos on real arguments

sin and cos on real arguments have a period of 2p1 .  When y is sufficiently large that y+2p1 has the same representation as y , then 2p1|!.0 x is 0 for all x>:y and sin x and cos x are meaningless. (sin x is 0 and cos x is 1 for all such x .)

The following demonstrates that y=:2^56 is the smallest number such that y+2p1 has the same representation as y .

```   ] y=: 2^56
7.20576e16
y+2p1
7.20576e16
y - y+2p1
0

hex y+0 2p1
4370000000000000
4370000000000000
hex (unhex '436',13\$'f')+0 2p1
436fffffffffffff
4370000000000000
```

2p1 |!.0 x is 0 for all x>:y :

```   ] x=. unhex ('437',12\$'0'),"1 0 '0123456789abcdef'
7.20576e16 7.20576e16 7.20576e16 7.20576e16 7.20576e16 7.20576e16 ...
hex x
4370000000000000
4370000000000001
4370000000000002
4370000000000003
4370000000000004
4370000000000005
4370000000000006
4370000000000007
4370000000000008
4370000000000009
437000000000000a
437000000000000b
437000000000000c
437000000000000d
437000000000000e
437000000000000f
x - {. x
0 16 32 48 64 80 96 112 128 144 160 176 192 208 224 240
2p1 |!.0 x
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
```

The J interpreter stops well short of the 2^56 limit. A limit error is signalled when 2p1 |!.0 x loses half of the number of bits of precision, that is, when the magnitude of x is greater than or equal to 2p1 * 2^53%2 where 53 is the number of bits in the significant (including the implied leading 1 bit). Thus:

```   thmax=: 2p1 * 2^53%2
thmax
5.96314e8

sin 6e8
|limit error: sin
|       sin 1000000000
cos _2e9
|limit error: cos
|       cos _2000000000
sin _
|limit error: sin
|       sin _
```

## Exponential

Mathematically, ^x is positive for all x . However, if ^x is smaller than the smallest positive IEEE number e , that is, if x < ^.e , then ^x is represented as 0 .

```   e
4.94066e_324
emin=: ^. e
emin
_744.44
_1074 * ^. 2
_744.44

^ emin
4.94066e_324

hex emin
c0874385446d71c3
x=. unhex '4' ,~ _1}. ,hex emin  NB. the next smaller number
hex x
c0874385446d71c4
^ x
0
```

At the other extreme, if ^x is larger than the largest IEEE number E , that is, if x > ^.E , then ^x is represented as _ (infinity).

```   E
1.79769e308
emax=: ^.E
emax
709.783
1024 * ^.2
709.783
hex emax,1024*^.2
40862e42fefa39ef
40862e42fefa39ef

^ emax
1.79769e308

hex emax
40862e42fefa39ef
x=. unhex 'f0' ,~ _2}. ,hex emax  NB. the next larger number
hex x
40862e42fefa39f0
^ x
_
```

Now we consider complex arguments. From the identity

```(^ x j. y) = (^x) * (cos y) j. (sin y)
```

If x<emin , ^x is (represented as) 0 and the product is 0 . Otherwise, if thmax<:|y , cos y and sin y are not accurate and a limit error is signalled.

If thmax>|y and x>:emin, the above identity is equivalent to:

```(^x j. y) = ((*cos y)*^x+^.|cos y) j. ((*sin y)*^x+^.|sin y)
```

The expression (*cos y)^x+^.|cos y is computationally more useful than (cos y)*^x for x>emax (similarly for the imaginary part).

```   ^ 716j1.57
7.17696e307j_
emax
709.783
^ 716
_

^ 716 + ^. cos 1.57
7.17696e307
^ 716 + ^. sin 1.57
_
```

The procedure can be modelled as follows:

```expcheck=: 3 : 0
'x y'=. +. y
if. x<emin do. 0j0 return. end.
assert. thmax>|y
if. x<:emax do. (^x) * (cos y) j. (sin y)
else. ((*cos y)*^x+^.|cos y) j. ((*sin y)*^x+^.|sin y) end.
)

expc 3j4
_13.1288j_15.2008
^ 3j4
_13.1288j_15.2008

expc 716j1.57
7.17696e307j_
^ 716j1.57
7.17696e307j_

expc 1j1e9
|assertion failure: exp
|   thmax>|y
^ 1j1e9
|limit error
|       ^1j1e9

expc _800j1e9
0
^ _800j1e9
0
```

## sinh and cosh

By definition,

```(sinh y) = 2 %~ (^y) - ^-y
(cosh y) = 2 %~ (^y) + ^-y
```

For large magnitude y, the computation comes down to 2%~^|y . From the Exponential section above it is seen that if emax<(|y)-^.2 , ^(|y)-^.2 exceeds the largest IEEE number E . That is, if (|y)>emax+^.2 , then sinh y and cosh y are _ .

```   emax2=: emax+^.2
emax2
710.476
E +&^. 2
710.476

sinh _1 1*emax2
_1.79769e308 1.79769e308
cosh _1 1*emax2
1.79769e308 1.79769e308

hex emax2
408633ce8fb9f87d
x=. unhex 'e',~}:,hex emax2  NB. the next larger number
hex x
408633ce8fb9f87e
sinh _1 1 * x
__ _
cosh _1 1 * x
_ _
```

For complex arguments, the identities

```(sinh z) = sin&.j. z    NB. A&S 4.5.7
(cosh z) = cos& j. z    NB. A&S 4.5.8
```

reduces sinh and cosh to the trigonometric versions. (See below.)

## sin and cos

```(sin x j. y) = ((sin x) * cosh y) j.   (cos x) * sinh y    NB. A&S 4.3.55
(cos x j. y) = ((cos x) * cosh y) j. - (sin x) * sinh y    NB. A&S 4.3.56
```

Since cosh y and sinh y can not be both 0 , for large magnitude x either the real part or the imaginary part (or both) of  sin x j. y and  cos x j. y are indeterminate, and a limit error should be signaled. Otherwise, the identities above apply.

The procedures can be modeled as follows:

```sincheck=: 3 : 0
'x y'=. +. y
assert. thmax>|x
((sin x)*cosh y) j. ( cos x)*sinh y
)

coscheck=: 3 : 0
'x y'=. +. y
assert. thmax>|x
((cos x)*cosh y) j. (-sin x)*sinh y
)

sincheck 3j4
3.85374j_27.0168
sin 3j4
3.85374j_27.0168

sincheck 2e9j_1
|assertion failure: sincheck
|   thmax>|x
sin 2e9j_1
|limit error: sin
|       sin 2e9j_1

sincheck 5j1000
__j_
sin 5j1000
__j_

sincheck 0j900
0j_
sin 0j900
0j_

coscheck 0j900
_
cos 0j900
_
```

## Collected Definitions

```sin  =: 1&o.
cos  =: 2&o.
sinh =: 5&o.
cosh =: 6&o.

hex  =: (4+#@\$) }. 2&(3!:3)
unhex=: (3!:2) @ ((}: 2 (3!:3) 0.5)&,) " 1

e    =: 2^_1074           NB. smallest positive finite number
E    =: (2-2^_52)*2^1023  NB. largest  positive finite number

thmax=: 2p1 * 2^53%2      NB. upper bound on arguments to sin and cos
emin =: ^.e               NB. ^y is 0 if y<emin;     _1074*^.2
emax =: ^.E               NB. ^y is _ if y>emax;      1024*^.2
emax2=: E +&^. 2          NB. cosh y is _ if y>emax2; 1025*^.2

expcheck=: 3 : 0
'x y'=. +. y
if. x<emin do. 0j0 return. end.
assert. thmax>|y
if. x<:emax do. (^x) * (cos y) j. (sin y)
else. ((*cos y)*^x+^.|cos y) j. ((*sin y)*^x+^.|sin y) end.
)

sincheck=: 3 : 0
'x y'=. +. y
assert. thmax>|x
((sin x)*cosh y) j. ( cos x)*sinh y
)

coscheck=: 3 : 0
'x y'=. +. y
assert. thmax>|x
((cos x)*cosh y) j. (-sin x)*sinh y
)
```