Sinusoids and harmonic oscillation

Harmonic oscillation is characterized by quantities which vary in such a way that each is proportional to its second derivative, the constants of proportionality all being equal and negative. Any solution to the differential equation used to encode this is expressible in terms of sinusoids which: rise smoothly to a peak, turn downwards thereafter, vary most rapidly while passing through zero, descend smoothly through a trough, and ascend once more to resume the cycle and repeat itself indefinitely, both backwards and forwards; and do so smoothly and symmetrically, as constrained by the differential equation.

Geometry gives us two standard sinusoids, Sin and Cos, which map angles to real scalars. I'll defer their formal definition to later: but, crucially, each is harmonic and, between them, they suffice to describe all harmonic oscillators. Note, however, that I deal with angles as having various units of measurement – notably the turn (360 degrees) and the radian (a bit over 57° or 1/(2.π) turns). The geometric and differential tool-kits each give us these sinusoids: on this page I intend to explore the implications of that (with at least half an eye on the issue of whether angles are dimensionless, or what the units of angle mean).

The Sinusoids

Sin and Cos are formally defined (along with Tan, Sec, coSec, coTan, which can all be inferred from Sin and Cos) in terms of ratios of sides of right-angle triangles, associating the ratio with one of the non-right angles. (Treating each side of a triangle as synonymous with its length …) If two sides, x and y, of a triangle meet in a right angle, the third, h satisfies h.h = x.x+y.y; if the angle between x and h is a and that between y and h is c, then Sin(a) = y/h while x/h = Cos(a). Since the sum of the angles in a triangle comes to a half-turn, and the angle between x and y is a quarter turn (a.k.a. right angle), the sum of a and c must be the remaining quarter turn; so if each is positive, each must be less than a quarter turn.

Symmetry of the description given implies Cos(c) = y/h and Sin(c) = x/h, implying Cos(c) = Sin(a) whenever a+c is a quarter turn. Thus far I've given meaning to Sin and Cos on angles (strictly) between 0.turn and turn/4. Take the corner at which x met h and use it as the origin for Cartesian co-ordinates: use the direction from this corner, along x, to the right angle as primary direction; as secondary, use the direction from there, along y, to where y meets h; through any point in the plane we may construct a straight line parallel to y; this must meet the continuation of x somewhere, making a right angle with it; with the ray from the origin to the given point we thus make a right-angle triangle with the ray as hypotenuse; the x co-ordinate of the point is the length of this triangle's side parallel to x, counted positive if this side is a segment of x or has x as a segment, otherwise negative; the y co-ordinate is the length of the side parallel to y, counted positive if the point is on the same side of the continuation of x as y, negative otherwise.

For any position in the positive quadrant, the triangle just given will yield the same co-ordinates as our original x, y, h, a, b triangle; so consideration of this one example suffices to describe all positions in the positive quadrant, each taken to be the c-corner. The standard representation of this corner in our Cartesian co-ordinates is [x,y].

Re-arranging the definitions of Sin(a) and Cos(a), we can write this representation as [h.Cos(a), h.Sin(a)] or equally as h.[Cos(a),Sin(a)]. Every position in the (interior – i.e. not on the axes – of the) positive quadrant can be expressed as h.[Cos(a),Sin(a)]: the ray from the origin to the position has length h and meets our original side, x, in an angle a. For positions outside the positive quadrant, we can construct the corresponding ray; it still has a length; we can divide each component by the length to get [c,s] for which c.c+s.s = 1, presenting c and s as values for Cos and Sin to take when given the angle between x and the ray.

Now the angle between x and a ray from the origin is synonymous with the rotation about the origin through that angle, which lines x up with the ray (and not, for instance, opposite to it). Either way, an angle has a magnitude – how large a fraction of the circle does it embrace – and a direction – the sense of rotation, clockwise or anti-clockwise. One can represent clockwise rotations as negative anti-clockwise ones, or indeed as (genuinely) anti-clockwise rotations through a turn minus the clockwise angle. Since rotation through an arbitrary whole number of turns (in either direction) is the identity, angles which (for a given choice of sense) differ by exact multiples of a turn are equivalent.

The addition on angles is effectively defined by composition of rotations; the composite is also a rotation, through an angle thus construed as the sum of the angles of the rotations composed. A rotation about the origin can be entirely described in terms of the position to which it maps any other point, so take [1,0] (i.e. our edge x, suitably scaled) and observe where a rotation takes it: the construction above applies directly, with h = 1, to give [x,y] = [Cos(a),Sin(a)] as the image of [1,0] under a rotation through a, at least when a is between turn and turn/4. From the structure of rotations, we can infer that [0,1] is then mapped to [−Sin(a),Cos(a)]. This specifies the rotation on two linearly independent directions in our (presumed) two-dimensional (real) space, hence determines it everywhere.

When we compose two rotations, through a and b, we map [1,0], via [Cos(b),Sin(b)] = Cos(b).[1,0] +Sin(b).[0,1], to Cos(b).[Cos(a),Sin(a)] +Sin(b).[−Sin(a),Cos(a)] = [Cos(b).Cos(a)−Sin(b).Sin(a), Cos(b).Sin(a)+Sin(b).Cos(a)] but the composite is just rotation through a+b, so this must be [Cos(a+b),Sin(a+b)], at least when a, b and a+b are all between zero (turns) and a quarter turn. Thus

Now write d for a+b, so b = d−a, and combine these to infer

Substituting various angles for b or d (e.g. c, for which a+c = turn/4) and even for a, we can obtain:

Thus the addition formula gives us the subtraction formula and, given the Pythagorean contribution, enables us to infer all the properties we knew before. When we look at the formulae for Cos and Sin of 3.a, we find that Cos(3.a) is zero whenever Cos(a) is zero and, otherwise, exactly when Sin(a) is ±1/2 (i.e. either square root of 1/4; implying Cos(a) = ±(√3)/2); likewise, Sin(3.a) is zero whenever Sin(a) is zero or Cos(a) is ±1/2. Now Sin(3.a) = 0 means 3.a is a whole number of half turns, aka an even number of quarter turns, and Cos(3.a) = 0 means 3.a is an odd number of quarter turns. These make a an even or odd number of twelfths of a turn, respectively, – so use a twelfth of a turn (a.k.a. 30 degrees) as unit for the present, calling it an hour; 3.hour = turn/4 – and obtain:

In either case, if the odd or even number is a multiple of three, a is an even or odd multiple of a quarter turn, 3.hour, hence the first conditions in the two cases, Sin(a) = 0 and Cos(a) = 0. Otherwise, the number isn't a multiple of 3, this first condition doesn't hold true and Cos(a) or Sin(a), as appropriate, must be ±1/2. This, combined with our knowledge that both Cos and Sin are positive between 0.hour and turn/4 = 3.hour, tells us Sin(1.hour) = 1/2 = Cos(2.hour). 3.hour is a quarter turn so Cos(3.hour+a) = −Sin(a) and Sin(3.hour+a) = Cos(a). Since 3+odd is even and 3+even is odd, we can infer:

Now, Pythagoras' contribution to that, h.h = x.x+y.y implies 1 = (x/h).(x/h) +(y/h).(y/h) = Cos(a).Cos(a) +Sin(a).Sin(a), i.e. (: Sin(a).Sin(a) +Cos(a).Cos(a) ←a :) is the constant function with value 1. Differentiating this, we obtain:

implying Sin'(x) = k(x).Cos(x), Cos'(x) = −k(x).Sin(x) for some ({linear ({reals}: :{angles})}: k :{angles}). For a given angle a and any mapping ({angles}:f:{reals}), in any neighbourhood in which f is differentiable, we obtain:

(: Sin(f(x)+a) ←x :)'
= (: f'(x).Sin'(f(x)+a) ←x :)
= (: f'(x).k(f(x)+a).Cos(f(x)+a) ←x :)

but, equally,

(: Sin(f(x)+a) ←x :)'
= (: Sin(f(x)+a) = Sin(f(x)).Cos(a) +Cos(f(x)).Sin(a) ←x :)'
= ( Cos(a).(Sin&on;f) +Sin(a).(Cos&on;f) )'
= Cos(a).(Sin&on;f)' +Sin(a).(Cos&on;f)'
= f'.(k&on;f).( Cos(a).(Cos&on;f) −Sin(a).(Sin&on;f) )
= (: f'(x).k(f(x)).Cos(f(x)+a) ←x :)

whence we infer k(f(x)+a) = k(f(x)), at least when f'(x).Cos(f(x)+a) isn't zero; and when it's zero a similar examination of (: Cos(f(x)+a) ←x :)' will yield k(f(x)+a) = k(f(x)) at least when f'(x).Sin(f(x)+a) isn't zero. Thus, for any angle a, any mapping ({angles}: f :{reals}), in any neighbourhood on which f is differentiable with non-zero derivative, k(f(x)+a) = k(f(x)).

For any angles b, d we can take f = (: x.d ←x :), a = b−d to show k(b) = k(f(1)+a) = k(f(1)) = k(d), since f is indeed differentiable: f(x) = x.d = x.f(1), f(x+y) = (x+y).f(1) = x.f(1) +y.f(1) = f(x) +f(y) so f is linear; hence equal to its own derivative (with respect to x); so f is differentiable, with constant derivative, everywhere. Note that f(x) = x.d, d = f(1) establishes a natural isomorphism between any additive domain and the collection of linear mappings from its scalings to the additive domain: (: x.d ←x :) is then synonymous with d.

Infer that Sin' = k.Cos and Cos' = −k.Sin with k constant. Thus Sin'' = −k.k.Sin and Cos'' = −k.k.Cos. For angles near 0.turn we can infer that, give or take terms tiny by comparison to a.a.a, Sin(a) = a.Sin'(0) + a.a.Sin''(0)/2 +a.a.a.Sin'''(0)/6 = a.k.Cos(0) −a.a.k.k.Sin(0)/2 −a.k.a.k.a.k.Cos(0)/6 = a.k −power(3, a.k)/6. Thus Sin(a) is slightly less than a.k, the fractional error being of order the square of a.k, so negligible in so far as a.k is small.

The position [Cos(a),Sin(a)] is close to [1,0] and both lie on the unit circle, so the straight line between these two positions is very close to the circle's arc between the same positions. The arc's length is a/turn times the circumference, 2.π, of our (unit) circle; let radian = turn/(2.π) to make the arc's length a/radian. We thus find a/radian and a.k are approximately equal to Sin(a) for all small angles, a, each approximation getting better the smaller a gets: we infer k = 1/radian, whence:

So differentiation measures angles in radians, even though the geometric symmetries of the sinusoids measure angles in turns.

This is the root cause of the widespread choice to use the radian as the definitive unit of angle: thus the functions provided by programming languages are sin = (: Sin(x.radian) ←x :{scalars}) and cos = (: Cos(x.radians) ←x :{scalars}). Now sin' = (: Sin'(x.radian).radian ←x :) = cos and cos' = −sin without the scalings: in terms of differential structure, the radian is the natural unit.

It remains that the turn is at least as natural – it encodes the periodicity and integer fractions of it appear in the symmetries of the sinusoids. Yet this requires factors of pi to describe in terms of radians, i.e. in terms of sin and cos.

The differential structure of sin and cos lead to the convenient formula exp(i.x) = cos(x) + i.sin(x) which arises when we combine an imaginary square root, i, of −1 with a function exp for which exp' = exp or, to be a little more pedantic, derivative(exp) = times&on;exp, where times = (: (: a*b ←b :) ←a :). The formula holds true when exp is further constrained to satisfy exp(0) = 1, yielding the exponentiation function commonly known as exp.

The structural mappings, multiplication and addition, on scalars intertwine with the exp mapping to give exp(a+b) = exp(a).exp(b) – without the constraint exp(0) = 1, the left side would say exp(0).exp(a+b), or the whole could be re-phrased exp(a).exp(b) = exp(c).exp(d) iff a+b = c+d. (The latter has some natural convenience for my preference for not presuming that 0 is available.)

This (taking i as name for a square root of −1) raises the question do pi and i have units – and, if so, in what contexts ?

Differentiating (: cos(x) + i.sin(x) ←x :) we get (: −sin(x) +i.cos(x) ←x :) = i.(: cos(x) +i.sin(x) ←x :) and differentiating (: exp(i.x) ←x :) we get (: i.exp'(i.x) ←x :) = i.(: exp(i.x) ←x :), so the two sides of the formula satisfy the same first order differential equation, f'=i.f, and agree on f(0) = 1.

So radians describe the differential structure of the sinusoids; yet turns describe their periodic structure. The place where the two meet is via exp when defined on the complex numbers, where it has periodicity 2.pi.i (as should be apparent from the formula, as 2.pi.radian = 1.turn = 360.degree is the period of Sin and Cos). This suggests to me that i may play a role in the story of units of angle; as may the differential operator. The reason I'm taking so much time and effort over this seemingly trivial matter is that angular momentum is deeply entwined – Dirac's constant, times a square root of −1, is the constant of proportionality between the differential operator (applied to a (spatially-described) wave-function) and the momentum operator (applied to the quantum state encoded by that wave-function). [Why do I keep mis-typing wave as weave ? ;^] So even if angles are dimensionless, there's some imaginary character (invariably entangled with a factor of 2.pi) to them which warrants handling in terms similar to how we handle dimensions – i.e. units.

I'm inclined to define turnExp = (: exp(2.pi.x) ←x :), which has derivative 2.pi.turnExp and period i, obtained from turnExp(x) = cos(x.turn) +i.sin(x.turn). In some important sense, sin and cos are describing a property of a linear space of real dimension 2 which is properly encoded by the complex numbers – by expressing trigonometry in terms of rotations, enlargements and one reflection. Their simple harmonic character (sin'=cos, cos'=−sin so, in particular, each is a solution of f''=−f, the defining equation of the simple harmonic oscillator) is encoded by a square root of −1 and a solution to f'=f. Given f' proportional to f with f(0)=1, we get f = exp if the constant of proportionality is 1; but if, instead, we constrain f to have period i we get f = turnExp, making the constant of proportionality 2.pi.

What's happening above is a choice, among the functions f for which f'/f is constant, of a canonical function: all others can be expressed in terms of any given one. The constraint f(0) = 1 arises naturally; then f(t) is f(1) raised to the power t, at least for integer t. The radian choice selects the candidate for which f'/f is also 1 and calls it exp; any other answer with f(0) = 1 will be (: exp(t.k) ← t :) for some complex k = f'/f. We can also ask for the candidate for which f(i) = i. Now, i is exp(i.π/2), so this selects f = (: exp(t.π/2) ← t :). [Exercise: what are the fixed points of exp and of other f with f'/f constant ?]

constraintk = f'/ff(1) = exp(k)periodfixed point
f'/f = 112.718 and a bit2.π.i0.318 +1.337i and a bit
f(−1) = −1; or
period = 2
π.i−12−1
f(i) = i; or
f(−i) = −i
π/24.81 and a little bit4.ii and −i
f(i) = −1π23.14 and a little bit2.i−0.1594 +0.5847i and a bit
f(i) = 12.π535 and nearly a halfi0.035 +1.2455i and a bit
period = 1; or
f(1) = 1
2.π.i110.15 +0.2133i and a bit
f'/f is a fixed point of f0.614 +0.681i and a bit 1.436 +1.16i and a bit5.0865 +4.588i and a bitf'/f
period = f'/f√(2.π.i) = (1+i).√π
≈ (1+i)×1.77
−1.17878 +5.766i and a bit(1+i).√π0.1796 +0.518i and a bit

I rather like the specifications in terms of fixed point, albeit they need not specify unique solutions (likewise the last constraint also has f'/f = −(1+i).√π as a solution). The f(1) = 1 case also has the constant solution f = (: 1←x :), among others (which is why I only list it as a variant on period = 1; which gives f(i) = exp(−2.π) = 0.001867 and a bit, while 1/f(i) = f(−i) = exp(2.π) = 535 and nearly a half). Asking for i = f(i) or −1 = f(−1) are among the simplest fixed point constraints available; the last two constraints are comparatively self-contained.

Python code for numerical solution of the exercise (with f(0) = 1):

from cmath import *

def fix(k, tiny=7e−17):
    """Returns a fixed point of (: exp(k.x) ← x :)

    i.e. some x with exp(k.x) = x.  Required first argument is k; optional
    second argument is a tolerance to within which you want the equation
    satisfied.  Default tolerance is just slightly larger than I found
        (exp(pi*.5j) −1j).real
    to be: it 'should be' zero but doesn't quite work out so.\n"""

    if k: x = 1j * k / abs(k) # a pure phase sideways to k's phase.
    else: return 1 # the unique fixed point of (: 1 = exp(0.x) ← x :)

    while True: # i.e. loop forever
        q = exp(k * x)
        d = q − x
        if abs(d) < tiny: return x # our only exit from the loop.
        x = x − d / (k * q − 1)

# Note that different initialization *could* produce different solutions.
# (Solving for f'/f a fixed point of f, i.e. k = fix(k), requires more effort …)

Of course, the nicest thing about orthodox exp is its power series:

wherein

which makes it computationally tractable; and involves no transcendental numbers in the computation (unlike the answers). All the other candidates above can, conveniently, be expressed in terms of it; each, notably, involves non-trivial scaling of the input. None the less, f(−1) = −1 provides the rather nice property f(even) = 1, f(odd) = −1 for even and odd integers; it is tantamount to using the half turn as unit of angle.


Valid CSS ? Valid HTML ? Written by Eddy.