APPENDICES

Equations: Some Basic Terms

There are two types of analytic equiations: explicit
and implicit.

Explicit equations are of the form y = f(x)
whereas implicit equations are of the form f(x,y)
= 0. The former is good for generating points while the
latter is good for testing to see if points satisfy the equation.
The drawback of the explicit form is that it is dependent on the
choice of coordinate axes and is ambiguous if there is more than
one y for a given x (such as y = sqrt(x)).

An alternative to the analytic form of equations is called
the parametric form, such as x = f(t), y =
g(t). This is good for generating points by varying the
parameter t.

Equations can be classified according to the terms contained
in them. Equations that only contain variables raised to a power
are polynomial equations. If the highest power
is one, then the equation is linear. If the
highest power is two, then the equation is quadratic.
If the highest power is three then it is cubic.
If the equation contains sins, cosines, log or a variety of other
functions, then it is called transcendental.

Continuity refers to how well behaved the
curve is in a mathematical sense. If, for a value arbitrarily
close to a value x0 the function is arbitrarily close to f(x0),
then it has positional, or zeroth order continuity
at that point. If the slope of the curve (or the first derivative
of the function) is continuous, then the function has first
order continuity. This is extended to all of the
functions derivatives.

If a curve is pieced together from individual curve segments
then one can speak of piecewise continuity if
each individual segment is continuous and the values at the
juntions of the curve segments are continuous.

If some data points are used to construct a curve the curve
can either pass through the points, in which case the curve is interpolating,
or the points can just be used to control the general shape of
the curve and the curve doesn't actually pass through the points,
in which case the curve is approximating.

An example: linear interpolation

Simple linear interpolation is given by

P(t) = (1-u)*P0 + u*P1

Notice that the interpolants, '1-u' and 'u', sum to one. This
property insures that the interplating curve (in this case a
straight line) falls within the convex hull of the geometric
entities being interpolated (in this case the convex hull is the
straight line itself). This is called the convex hull
property.

Using more general notation, the equation above can be
rewritten as:

P(u) = F0(u)*P0 + F1(u)*P1

When written in this form, F0 and F1 are referred to as blending
functions. The form is referred to as the geometric
formulation because the geometric information, in this
case P0 and P1, are explicit in the equation.

The linear interpolation equation can be rewritten as:

P(u) = (P1-P0)*u + P0

Again, using more general notation, the equation is of the
form: P(u) = a1*u+ a0 in which the terms are collected into
coefficients of the variable to some power. This is the general
polynomial form referred to as the algebraic formulation.

Alternatively, both of these forms can be put in a matrix
representation . The geometric formulation becomes:

T T
P(t) = [ F0(u) F1(u) ] [P0 P1] = F B

T
(B denotes the transpose of matrix B.)

The algebraic form becomes:

T T
P(t) = [u 1] [a0 a1] = U A

In the fully expanded form, we have:

T T T T
P(u) = [u 1] [ -1 1 ] [P0 P1] = U M B = F B = U A
[ 1 0 ]

Most all of the forms for equations that we will deal with
can be written in the format above. Depending on the actual curve
type, the U (variable matrix), M (coefficient matrix) and P
(geometric information matrix) matrices will have different
information in them.

Reparameterization by ArcLength

It should be noted that, in general, there is not a linear
relationship between changes in the parameter u and the distance
traveled along a curve (its arclength). It does happen to be true
in the example above concerning a straight line and the parameter
u. However, as Mortenson
points out, there are other equations which trace out a straight
line in space which are fairly convoluted in their relationship
between changes in the parameter and distance traveled. The
non-linear relationship is almost always true in any
parameterized curve unless special care is taken to ensure
linearity. (See
Reparameterization by Arc Length in the Chapter on
Techniques to Aid Motion Specification)

Hermite Interpolation

Hermite formulation produces a cubic polynomial through the
points. The user needs to supply beginning/end tangents as well
as end points.

T
P(u) = U M B

where

3 2
U = [u u u 1]

M = matrix of coefficients found by solving the cubic equation
with constraints

P = [P0 P1 P0' P1']

where P0, P0' are the beginning point and tangent vector,
respectively, and P1, P1' are the ending point and tangent
vector, respectively.

Continuity between beginning and ending tangent vectors of
connected segments is ensured by merely using the ending tangent
vector of one segment as the beginning tangent vector of the next
one. For a single cubic Hermitian curve, we have:

If trying to put a Hermite curve through a large number of
points, requiring the user to specify all of the needed tangent
vectors can be a burden. There are several techniques to get
around this. One such technique is to enforce second degree
continuity. This requirement provides enough constraints so that
the user doesn't have to provide interior tangent vectors; they
can be calculated automatically. See Rogers
and Adams' book (p.119-123 in the first edition), or see Mortenson's book (p. 98-112) for
an alternative formulation.

Catmull-Romm Spline

The Catmull-Romm curve can be viewed as a Hermite curve in
which the tangents at the interior control points are
automatically generated according to a relatively simple
geometric procedure (as opposed to the more involved numeric
techniques referred to above). For each interior point, p[i], the
tangent at that point, p'[i], is computed as:

p'[i] = (1/3)*(p[i+1]-p[i-1])

Catmull-Romm spline is a specific type of the cardinal
splines. In cardinal splines, the factor of 1/3 is user
specified. In addition, the beginning and end tangent vectors can
also be automatically determined by various simple geometric
means. For example,

p'[0] = (1/3)(p[1]-p[0])
or
p'[0] = (1/2)*(2*p[1]-p[2]-p[0])

and these definitions for the tangents can be put in the B
matrix of the equation for Hermite interpolation given above.

Blended Parabolas

A cubic curve can also be created in the middle section of
four control points in the following way. Take the first three
points, P1, P2, P3, and fit a parabola, P(u), through them using
the following contraints: P(0) = P1, P(1/2) = P2, P(1) = P3. Take
the last three pints, P2, P3, P4, and fit a parabola, R(u),
throught them using similar constraints: R(0) = P2, R(1/2) = P3,
R(1) = P4. Between points P2 and P3 the two parabolas overlap.
Reparametrize this region into the range [0,1] and linearly blend
the two parabolic segments:

For interpolating a list of points, interior segments are
calculated using the equation above. End conditions can be
effected by making parabolic arcs at the very beginning and very
end.

This form assumes that all points are equally spaced in
parametric space. Often it is the case that even spacing is not
present. In such cases, relative cord length can be used to
estimate parametric values.

Bezier Interpolation/Approximation

The are several ways to look at Bezier curves. One way to look
at cubic Bezier curves is to say that it is similar to the
Hermite formulation, except instead of using tangent vectors, the
Bezier formulation uses auxiliary control points to define
tangent vectors. A cubic curve is defined by four points: P0, P1,
P2, P3. P0 and P3 are the beginning and end points, respectively,
of the curve. P1 and P2 are the interior control points which
define the beginning and ending tangents.

Continuity between adjacent Bezier segments can be controlled
by colinearity of the control points on either side of the shared
beginning/end point of the two curve segments where they join.

In addition, the Bezier curve formulation allows one to
define a curve of arbitrary order. If three interior control
points are used, then the resulting curve will be quartic; if
four interior control points are used, then the resulting curve
will be quintic. See Mortenson
for a discussion. For a single cubic Beziercurve:

Bspline Approximation/NURBS

B-splines are a more general formulation which includes Bezier
as a special case. The formulation decouples the number of
control points from the degree of the resulting polynomial with
the knot vector. Values in the knot vector also allow
specification of control points to be interpolated or just
approximated by the curve. For a single cubic B-spline curve:

NURBS, which we won't get into here, are even more flexible
than basic B-splines. Nonrational Bsplines and Bezier curves are
special cases of NURBS. NURBS allow for exact representation of
circular arcs whereas Bezier and nonrational Bsplines don't.

Bias and Tension

Considered by some to be more intuitive control of the curves.
P(t) = TMP

Physically-based motion is really just a limited simulation of
physical reality. This can be as simple or as complex as the
implementor desires. Below are some of the equations of
importance in simple physics simulation.

Position, Velocity, Acceleration

The fundamental equation relating position, distance, and
speed is

distance = speed * time

This can be used to control the positioning of an object for
any particular frame of the animation because this is tied
directly to time

time = frame-number * time-per-frame

The average velocity of a body is the distance moved divided
by the time it took to move, and is given by

v = (distance between s(t1) and s(t2))/(t2-t1) (B.1)

where s(t) is a function that gives the position of the
object at time t. Notice that the units of velocity is distance
per time, e.g., ft/sec.

The instantaneous velocity is determined by taking equation
(B.1) and moving t2 closer and closer to t1. In the limiting case
this becomes the derivative of the function s(t) with respect to
time.

Similarly, the average acceleration of an object is the
change in velocity divided by the time it took to effect the
change, and is given by

a = (v(t2)-v(t1))/(t2-t1)

where v(t) is a function that gives the velocity of the
object at time t. Notice that the units of acceleration is
velocity per time or distance per time per time, e.g., ft/sec2.

In the same way, instantaneous acceleration is the derivative
of v(t) with respect to t so that we now have

a(t) = v'(t) = s''(t)

For example, in the case of motion due to gravity, we have:

s(t) = (1/2)gt**2
v(t) = gt
a(t) = g

where g is the acceleration due to gravity, a constant that
has been measured to be 32 ft/sec2.

Circular Motion

Circular motion is an important type of motion in physics and
comes up when discussing cases from planets to robotics. A simple
way to describe it is to specify position based on polar
coordinates. The position of a particle orbiting the origin a
distance r from it can be described as:

p = (r*cos(theta))i + (r*sin(theta))j (B.2)

where i and j are orthonormal unit vectors (at right angles
to each other and unit length). p is the positional vector of the
particle. In a constant radius curcular orbit, the distance r is
constant, but it is important to keep in mind that p, being a
vector, is not.

Uniform circular motion means that theta changes at a
constant rate which we call omega:

d(theta)/dt = omega

or theta = omega*time

If theta increases with time, the quantity omega is positive
as is called the angular speed of the object. If theta is
measured in radians and time in seconds, then omega is measured
in radians per second.

The angular speed, omega, is related to the time T that it
takes to complete one revolution. T is called the period and is
related by the following equation:

2PI = omega*T

since omega moves through 2PI radians in one revolution.

Given (B.2), we can take the derivative of it with respect to
time to get the acceleration. After substituting (B.3) and
calling it p(t), using w to represent omega and taking the
derivative, we get:

v(t) = dp/dt = (-r*w*sin(w*t))i + (r*w*cos(w*t))j (B.4)

Notice that the velocity vector, v(t), is perpendicular to
the position vector p(t) (this can be demonstrated by taking the
dot product of the two vectors v(t) and p(t) to get zero. By
computing the length of v(t) we arrive at:

v = r*w

showing that the velocity is independent of t and is,
therefore, constant.

Notice, however, that a constant circular motion still gives
rise to an acceleration. We take the derivative of (B.4) to get:

a(t) = dv/dt = -w2*((r*cos(w*t))i + (r*sin(w*t))j)
= w**2*r

using v = r*w, we have

a(t) = v**2/r

which is called the centripetal acceleration which is
directed radially inward and has constant magnitude.

For any particle in a rigid mass undergoing a rotation, that
particle is undergoing the same rotation about its own center
and, in addition, because it may be displaced from the center of
rotation, it is undergoing an instantaneous positional
translation.

Newton's Laws of motion

First Law: If no force is acting on an object, then it will
not accelerate. It will maintain a constant velocity.

This is the principle of inertia.

Second Law: The change of motion of an object is proportional
to the forces applied to it.

F = m*a

The change in motion (acceleration) of an object due to a
force is proportional to the object's mass. That is, the change
in motion involved not only the body's velocity but also its
mass. Stated another way this law becomes

F = d(m*v)/dt

where the quantity mv is the object's momentum. In this form
the law states 'Force is equal to the rate of change of momentum.

It is important to note that the force F used here is
considered to be the sum of all external forces acting on an
object and that these equations really represent three sets of
equations, one for each coordinate:

Fx = m*ax
Fy = m*ay
Fz = m*az

Third Law: To every action there is always opposed an equal
and opposite reaction.

Inertia and inertial reference frames

An inertial frame is a frame of reference (e.g., coordinate
system) in which the principle of inertia holds. Any frame that
is not accelerating is an inertial frame. An observer in an
inertail frame deduces the same laws of motion and has no way of
determining whether he is at rest or moving in an absolute sense
(but then, what is moving in an absolute sense?)

torque

The tendency of a force to produce rotation is called torque.
In 2D, we have:

It is important to note that torque refers to a specific axis
of rotation. The same force can exert different torques about
different axes of rotation. The torque vector for a rigid body is
position independent. That is, for any particle in a rigid mass,
the particle is undergoing that same torque

Equilibrium: balancing forces

In the absense of acceleration, there is equilibrium. Study of
bodies in equilibrium is statics. In such as case the net force
and torque acting on a body vanishes:

SUM of F = 0
SUM of t = 0

There may be forces present, but the vector sum is equal to 0.

Gravity

Law of universal gravitation is

F = G*m1*m2/r**2 where G = 6.67x10-11 m3/kg sec2

The resultant force on a mass is the vector sum of the
individual forces.

When two spherically symmetric objects are not touching, each
acts as if all its mass were concentrated at its center.

The force of gravity from a mass distributed symmetrically
over a spherical shell is zero anywhere inside the shell.

Centripetal force

Centripetal force is any force which is directed inward,
toward the center of an object's motion (cetripetal means 'center
seeking'). In the case of a body in orbit, gravity is the
cetripetal force that holds the body in the orbit. Consider a
body, such as the moon, with mass Mm and a distance r away from
the earth. The unit vector from the earth to the moon is R. The
earth has a mass Me. The acceleration induced by a circular orbit
always points toward the center (centripetal) and has magnitude
v2/r:

a = (-v**2/r)R

The force inducing the acceleration is gravity and the force
of gravity is:

F = (-GMmMe/(r*r))R

and it causes the moon to have acceleration

a = F/Mm

or a = (-GMe/(r*r))R

using the original equation for acceleration of circular
motion and canceling terms we can solve for the velocity of the
moon which compensates for the effect of gravity and establishes
an orbit to be:

v = sqrt(GMe/r) Similar analyses can be carried out when
considering a ferris wheel, an airplane banking to make a turn or
the motion of a conical pendulum.

Contact Forces: Hook's Law

While gravity is one of the four fundamental forces (gravity,
electromagnetic, strong, and weak) which act at a distance. A
second important category is contact forces.The tension in a wire
or rope is an example of a contact force, as is the compression
in a rigid rod. These forces arise from the complex interactions
of electric forces which tend to keep atoms a certain distance
apart. The empirical law which governs such forces, and which is
familiar when dealing with springs, is Hook's Law:

F = -kx where x = change in lenght from the equilibrium value of the spring
k = spring constant
F = restoring force of sping to return to rest position

The constant k is a measure of the stiffness of the spring;
the larger k is, the more sensitive the spring is to motion away
from the rest position.

Another example of contact force is the result of repulsion
between any two objects pressed against each other, known as the
normal force. It arises from the repulsion of the atoms of the
two objects, is always perpendicular to the surfaces and its
magnitude is proportional two how hard pressed the two objects
are. Other important examples of contact forces are friction and
viscosity.

Friction

Frictional forces from a surface do not exceed an amount
proportional to the normal force exerted by the surface on the
object:

fs <= µsn where µs="coefficient" of static friction n="normal" force

µs varies with the two materials that are in contact. The frictional forces acting between surfaces at rest with respect to each other are called forces of static friction. For a force applied to an object sitting on another object that is parallel to the surfaces in contact, there will be a specific force at which the block starts to slip. As F increases from zero up until that threshold force, the lateral force is conteracted by an equal force of friction in the opposite direction. Once the object begins to move Kinetic friction acts on the object and approximately obeys the empirical law:

fk = µkN where µk = coefficient of kinetic friction N =
normal force

Kinetic friction is typically less than static friction. The
force of kinetic friction is always opposite to the velocity of
the object.

Motion in a resistive medium

The resistive force is viscosity. It is another contact force and extremely difficult to model accurately. When an object moves at low velocity through a fluid, the viscosity is approximately proportional to the velocity:

The coefficient of viscosity, n, decreases with increasing
temperature for liquids, and increases with temperature for
gases.

Stokes' law for a sphere of radius R gives:

K = 6¹R

An object dropping through a liquid attains a constant speed,
called the limiting or terminal velocity, at which gravity,
acting downward, and the viscous force, acting upward, balance
each other (no acceleration):

mg = 6¹Rnv
v = mv/6¹Rn

In a viscous medium, heavier bodies fall faster than lighter bodies. For sphereical objects falling a low velocity in a viscous medium, not necessarily at terminal velocity:

mdv/dt = mg - 6¹Rnv

centrifugal force

Consider a frame of reference that rotates in uniform circular motion with respect to an inertial frame, a post, because it is held at a constant distance by a rope. Relative to the inertial frame, each point in the uniformly rotating frame has centripetal acceleration:

a = -(v**2/r)R where
r is the distance the point is from the axis of rotation
R is the unit vector from the inertial frame to the rotating frame
v is the speed of the point

And the tension in the rope supplies the force necessary to produce the centripetal acceleration. Relative to the rotating frame (not an inertial frame), the frame itself does not move and therefore to counteract the force supplied by the rope is the centrifugal force:

Fc = -ma = mw**2rR

Work and Potential Energy

For a constant force of magnitude F moving an object a distance h parallel to the force, the work W done by the force:

W = Fh

If a mass m is lifted up so that it doesn't accelerate, then
the lifting force is equal to the weight (mg) of the object.
Since mg is constant, the work done to raise the object up to a
height h is

W = mgh

Energy that a body has by virtue of its location is called
potential energy. The work in this case is converted into
potential energy.

Kinetic Energy

Energy of motion is called kinetic energy. For a falling body that started at a height h,

v**2 = 2gh
(1/2)mv**2 = mgh
K = (1/2)mv**2

Conservation of energy

Power

Power is defined to be the rate at which work is done with respect to time.
P = dW/dt

Conservation of momentum

In a closed system, momentum (mass times velocity) is conserved:
d(m1v1 + m2v2 + ... + mnvn)/dt = 0
m1v1 + m2v2 + ... + mnvn = const.
The center of mass r of a system of particles moves under external forces as if all the mass were concentrated at r and all forces acted at r.

Impluse: Collision Forces and Times

[Need collision response material]

oscillatory motion

Stability of an object subject to a linear restoring force:

F = -kx

Associated with this force is the potential energy:

U = (1/2)kx2

These systems have the property that if they are distrubed
from equilibrium, the restoring force that acts on them tends to
move them back into equilibrium. When disturbed from equilibrium,
they tend to overshoot that point when they return, on account of
inertia. Then the restoring force acts in the opposite direction,
trying to return the system to equilibrium. The result is that
the system winds up oscillating back and forth like a mass on the
end of a spring.

We have the following:

F = ma
F = -kx
a = d2x/dt2
m d2x/dt2 = -kx
d2x/dt2 = -kx/m

is a differential equation satisfied by the displacement x(t).

Damping

We can model the damping force after Stokes' law, in which resistance is assumed to be proportional to velocity. This is usually valid for oscillations of sufficiently small amplitude. The damping force opposes the motion and can be written as:

Fd = -½*dx/dy

where ½ is aconstand called the damping coefficient. Adding
the damping force to the spring force, Newton's second law gives
us:

m*d2x/dt2 = -k*x - ½*dx/dt

After dividing through by m and rearranging:

d2x/dt2 + §*dx/dt + w20*x = 0

where § = ½/m and w20 = k/m.

If there is a sping force but no damping, the general
solution can be written as C*cos(w0*t+¯0). If there is damping
but no spring force, the general solution turns out to be x =
C*e-§t+D with C and D constant. If both the spring force and
damping force are present, the solution is more compicated.

Rotational Dynamics for Rigid Bodies

overview

This is a collection of techniques and pointers to other sources that should help the reader understand and implement numerical techniques that are useful for computer animation algorithms. This is intended just to give the reader a starting point.

There is a hierarchy of complexity of numerical problems. In
increasing order of complexity:

We will treat each of these in turn to give an overview of various approaches to handle problems in these domains. Usually the techniques will, by approximations, reduce a problem from one domain to a simpler domain until there is a procedure for determining a solution in the simpler domain.

Linear algebraic equations can be expressed in matrix form,
for example, the familiar

Ax = b

where A and b are a constant matrix and vector respectively
and x is a variable vector which is solved for.

Non-linear algebraic equations involve no derivatives but are
more complex than simple linear equations. For example,
polynomials fall into this class. Some non-linear equations can
be solved directly, for example, a quadratic polynomial can be
solved using the quadratic equation. Others can be converted into
linear algebraic equations by using simplifying assumptions and
then solved using the linear techniques.

Ordinary differential equations (ODEs) involve the derivative
of only one variable. Some ODEs can be solved directly. ODEs can
be converted into non-linear agebraic equations using simplifying
assumptions and solved from there.

Partial differential equations (PDEs) involve functions of
more than one variable and derivatives of more than one variable.
There are techniques which convert PDEs to ODEs.

Taylor series and differencing

A common useful equation that comes up is Taylor's series expansion:

f(x+¶) = f(x) + f'(x)¶ + (1/2)f''(x)¶2 + ...

This equation comes in handy when the current situation is
known or can be approximated, and future conditions must be
calculated. An example is the case where the motion of an object
is controlled by physical laws and f(x) is the position of that
object, f'(x) is the velocity and f''(x) is the acceleration
(e.g., under gravity). The position of the object at the next
instant of time (i.e., x+¶) is to be computed, acceleration is
given, for example, by gravity and velocity can be approximated
by adding the acceleration to the last velocity.

Linear equation solvers

Systems of linear equations arise, for example, when considering a mesh of mass points connected by springs for flexible body animation. Each pair of mass points considered to be connected by a spring gives rise to a linear equation of the form F = -kx.

Linear systems of equations can be written

Ax = b where A is a two-dimensional matrix of coefficients
x is a vector of unknowns
b is a vector of constants

The most obvious approach to this problem may appear to be to calculate the inverse of A and multiply both sides by it:

x = A-1b

The inverse of A, A-1, can be computed by Gauss-Jordan
elimination with either partial or full pivoting. It is more
efficient, however, to do Gaussian elimination until an upper
triangular matrix results and follow with backsubstitution. Both
of these solutions require that the right hand side of the
original equations be known in advance. A better approach is to
decompose A into an upper triangular matrix and a lower
triangular matrix:

A = LU where L is lower triangular U is upper triangular

so that

Ax = (LU)x = L(Ux) = b

Note that the LU decomposition of A can be done without
knowing the b vector. The reason for doing this is that solving a
system of equations with a triangular matrix is accomplished by
simple substitiution. So we first solve

Ly = b

called forward substitution, and then

Ux = y

called backsubstitution.

The only thing left to do is specify how the LU decomposition
is to be done. In order to keep the calculation stable, the
decomposition should be done with pivoting so that a permutation
of A is actually computed.

Singular value decomposition is a powerful set of techniques
for solving most linear least squares problems and for dealing
with systems of equations are either singular or else numerically
very close to being singular. The singular value decomposition of
a matrix which is MxN, where the number of rows M is greater than
or equal to the number of columns N, is an MxN column-orthogonal
matrix U, an NxN diagonal matrix W with positive or zero
elements, and the transpose of an NxN orthogonal matrix V.
Zeroing out small values in W essentially well-conditions the
resulting set of linear equations and very often produces a
better solution that both the direct-method solution and the SVD
solution where the small w's are left nonzero.

non-linear algebraic equations

Non-linear functions arise in a variety of places, most notably when dealing with higher order surfaces.

In a few cases, as mentioned earlier in this Appendix,
non-linear algebraic equations have analytic solutions, such as
the quadratic equation. However, most of the time analytic
solutions are not known for equations of this type. Simlar to the
linear case, a simple way to solve non-linear equation is by
functional iteration. If you have an equation in the form of

y = f(y)

you can use an initial guess and use it to generate the next
estimate to the answer:

yk+1 = g(yk)

However, functional iteration has poor convergence
properties.

A better technique for solving non-linear functions is
Newton's method which tries to find zero of non-linear function
by truncating the Taylor's series after two terms. !!!

There are other non-linear equation solvers, such as the
variable metric method. !!!

Such hill-climbing techniques are applicable in function of
two variables because, from a given point on a curve there are
only two directions you can proceed. With higher dimensional
functions, the problem is much more complicated because a surface
is defined and there are an infinite number of directions to
proceed from any one point.

ordinary differential equations

Problems involving ordinary differential equations (ODEs) can always be reduced to the study of sets of first order differental equations. For example, the secnd order equation

d2y/dx2 + q(x)dy/dx = r(x)

can be rewritten as two first order equations

dy/dx = z(x) dz/dx = r(x) - q(x)z(x)

where z is a new variable. The usual choice for new variables
is to let them be derivatives of eadch other and of the original
variable. Occassionally, it is useful to incorporate some other
factors into their definition fro the purpose of mitigating
singular behavior that could result in overflows or increased
roundoff errors.

The generic problem in ordinary differential equations is
thus reduced to the study of a set of N coupled first order
differential equations.

dyi(x)/dx = fi(x,y1,...,yN) i = 1,...,N

The nature of the problem's boundary conditions are crucial
in determining how to attack the problem. Boundary conditions are
algebraic conditions on the values of the functions yi. In
general they can be satisfied at discrete specified points, but
do not hold between those points, i.e., are not preserved
automatically by the differential equations. Boundary conditions
divide into two broad categories:

In initial value problems all the yi are given at some
starting value xs, and it is desired to find the yi's at
some final point xf, or at some discrete list of points
(for example, at tabulated intervals).

in two-point boundary value problems, on the other hand,
boundary conditions are specified at more than one x.
Typically, some of the conditions will be specified at xs
and the remainder at xf.

The underlying idea of any routine for solving the initial value problem is always to rewrite the dy's and dx's as finite steps Æy and Æx and multiply the equations by Æx. In the limit of making the stepsize very small, a good approximation to the underlying differential equation is achieved. Literal implementation of this procedure results in Euler's method whcih is not recommended for any practical use.

Three major types of practical numerical methods for solving
initial value prolems for ODEs are:

Runge-Kutta

Richardson extrapolation

predictor-corrector methods

partial differential equations

Partial differential equations are at the heart of many computer analyses or simulations of continuous physical systems, such as fluids, electromagnetic fields, the human body, tec. Partial differential equations are classified into three categories, hyperbolic, parabolic and elliptic, on the basis of their characteristics, or curves of information propagation. The prototypcial example of a hyperbolic equations is the one-dimensional wave equation:

¶2u/¶t2 = v2(¶2u/¶x2)

where v = constant is the velocity of wave propagation. The
prototypical parabolic equation is the diffusion equation:

¶u/¶t = -¶(D(¶u/¶x)/¶x

where D is the diffusion coefficient. The prototypical
elliptic equation is the Poisson equation:

¶2u/¶x2 + ¶2u/¶y2 = p(x,y)

where the source term p is given. If the source term is equal
to zero, the equation is Laplace's equation.

The first two are initial value problems which describe how
u(x,t) propagates itself forward in time. The last one is a
boundary value problem where the task is to find a single static
function u(x,y) which satisfies the equation within some (x,y)
region of interest and which has some desired behavior on the
boundary of the region of interest.

Components

Length

|A| = Ã(Ax2 + Ay2 + Az2)

Dot Product

A.B = |A||B|cos¿ where ¿ is the angle from A to B, 0 <= ¿ <="180" a.b="Ax*Bx" + ay*by + az*bz a.a="|A|2" a.b="0" if and only if a and b are perpendicular (i.e., cos(¿)="cos(90)" or a="0" or b="0" a.i="Ax" a.j="Ay" a.k="Az"

Cross Product
| i j k |
A x B = | Ax Ay Az |
| Bx By Bz |
A x B = (AyBz - AzBy)i + (AzBx - AxBz)j + (AxBy - AyBx)k
A x B = -(B x A)
the direction of A x B is perpendicular to both and is determined by the right hand rule (if A and B are in right-hand space)
|A x B| = |A||B|sin¿ where ¿ is the angle from A to B, 0 <= ¿ <="180" |a x b|="0" if and only if a and b are parallel (i.e., sin(¿)="sin(0)" or a="0" or b="0"