Clustering: Self-Organizing Map (SOM)

The Self-Organizing Map (SOM) is a clustering algorithm that is used to map a multi-dimensional dataset onto a (typically) two-dimensional surface. This surface (a map) is an ordered interpretation of the probability distribution of the available genes/samples of the input dataset. SOMs have been used extensively in many domains, including the exploratory data analysis of gene expression patterns.
There are two particularly useful purposes for this: visualization and cluster analysis. Visualization has typically been a difficult matter for high-dimensional data. SOMs can be used to explore the groupings and relations within such data by projecting the data on to a two-dimensional image that clearly indicates regions of similarity. Even if visualization is not the goal of applying SOM to a dataset, the clustering ability of the SOM is very useful.

MPEG 7 Standard

MPEG-7, formally named “Multimedia Content Description Interface”, is a standard for describing the multimedia content data that supports some degree of interpretation of the information’s meaning, which can be passed onto, or accessed by, a device or a computer code. MPEG-7 is not aimed at any one application in particular; rather, the elements that MPEG-7 standardizes support as broad a range of applications as possible.

The elements that MPEG-7 standardizes provide support to a broad range of applications (for example, multimedia digital libraries, broadcast media selection, multimedia editing, home entertainment devices, etc.). MPEG-7 will also make the web as searchable for multimedia content as it is searchable for text today. This would apply especially to large content archives, which are being made accessible to the public, as well as to multimedia catalogues enabling people to identify content for purchase.

OBJECTIVES OF MPEG-7 STANDARD

The MPEG-7 standard aims at providing standardized core technologies allowing description of audiovisual data content in multimedia environments. Audiovisual data content that has MPEG-7 data associated with it, may include: still pictures, graphics, 3D models, audio, speech, video, and composition information about how these elements are combined in a multimedia presentation (scenarios). Special cases of these general data types may include facial expressions and personal characteristics .

Geometric meaning of SVD

From Wikimedia

Because U and V are unitary, we know that the columns u₁,...,u_m of U yield an orthonormal basis of K^m and the columns v₁,...,v_n of V yield an orthonormal basis of Kⁿ (with respect to the standard scalar products on these spaces).

The linear transformation T :Kⁿ → K^m that takes a vector x to Mx has a particularly simple description with respect to these orthonormal bases: we have T(v_i) = σ_i u_i, for i = 1,...,min(m,n), where σ_i is the i-th diagonal entry of Σ, and T(v_i) = 0 for i > min(m,n).

The geometric content of the SVD theorem can thus be summarized as follows: for every linear map T :Kⁿ → K^m one can find orthonormal bases of Kⁿ and K^m such that T maps the i-th basis vector of Kⁿ to a non-negative multiple of the i-th basis vector of K^m, and sends the left-over basis vectors to zero. With respect to these bases, the map T is therefore represented by a diagonal matrix with non-negative real diagonal entries.

To get a more visual flavour of singular values and SVD decomposition —at least when working on real vector spaces— consider the sphere S of radius one in Rⁿ. The linear map T maps this sphere onto an ellipsoid in R^m. Non-zero singular values are simply the lengths of the semi-axes of this ellipsoid. Especially when n=m, and all the singular values are distinct and non-zero, the SVD decomposition of the linear map T can be easily analysed as a succession of three consecutive moves : consider the ellipsoid T(S) and specifically its axes ; then consider the directions in Rⁿ sent by T onto these axes. These directions happen to be mutually orthogonal. Apply first an isometry v* sending these directions to the coordinate axes of Rⁿ. On a second move, apply an endomorphism d diagonalized along the coordinate axes and stretching or shrinking in each direction, using the semi-axes lengths of T(S) as stretching coefficients. The composition d o v* then sends the unit-sphere onto an ellipsoid isometric to T(S). To define the third and last move u, just apply an isometry to this ellipsoid so as to carry it over T(S). As can be easily checked, the composition u o d o v* coincides with T.

Endomorphism

From Wikipedia

Endomorphism

In mathematics, an endomorphism is a morphism (or homomorphism) from a mathematical object to itself. For example, an endomorphism of a vector space V is a linear map ƒ: V → V and an endomorphism of a group G is a group homomorphism ƒ: G → G, etc. In general, we can talk about endomorphisms in any category. In the category of sets, endomorphisms are simply functions from a set S into itself.

In any category, the composition of any two endomorphisms of X is again an endomorphism of X. It follows that the set of all endomorphisms of X forms a monoid, denoted End(X) (or End_C(X) to emphasize the category C).

An invertible endomorphism of X is called an automorphism. The set of all automorphisms is a subgroup of End(X), called the automorphism group of X and denoted Aut(X). In the following diagram, the arrows denote implication:

automorphism		isomorphism
endomorphism		(homo)morphism

Any two endomorphisms of an abelian group A can be added together by the rule (ƒ + g)(a) = ƒ(a) + g(a). Under this addition, the endomorphisms of an abelian group form a ring (the endomorphism ring). For example, the set of endomorphisms of Zⁿ is the ring of all n × n matrices with integer entries. The endomorphisms of a vector space, module, ring, or algebra also form a ring, as do the endomorphisms of any object in a preadditive category. The endomorphisms of a nonabelian group generate an algebraic structure known as a nearring.

Operator theory

In any concrete category, especially for vector spaces, endomorphisms are maps from a set into itself, and may be interpreted as unary operators on that set, acting on the elements, and allowing to define the notion of orbits of elements, etc.

Depending on the additional structure defined for the category at hand (topology, metric, ...), such operators can have properties like continuity, boundedness, and so on. More details should be found in the article about operator theory.

Isometry

From Wikipedia

In mathematics, an isometry, isometric isomorphism or congruence mapping is a distance-preserving isomorphism between metric spaces. Geometric figures which can be related by an isometry are called congruent.

Isometries are often used in constructions where one space is embedded in another space. For instance, the completion of a metric space M involves an isometry from M into M', a quotient set of the space of Cauchy sequences on M. The original space M is thus isometrically isomorphic to a subspace of a complete metric space, and it is usually identified with this subspace. Other embedding constructions show that every metric space is isometrically isomorphic to a closed subset of some normed vector space and that every complete metric space is isometrically isomorphic to a closed subset of some Banach space.

Definitions

The notion of isometry comes in two main flavors: global isometry and a weaker notion path isometry or arcwise isometry. Both are often called just isometry and one should determine from context which one is intended.

Let X and Y be metric spaces with metrics d_Y and d_X. A map ƒ : X → Y is called distance preserving if for any x,y ∈ X one has

A distance preserving map is automatically injective.

A global isometry is a bijective distance preserving map. A path isometry or arcwise isometry is a map which preserves the lengths of curves (not necessarily bijective).

Two metric spaces X and Y are called isometric if there is an isometry from X to Y. The set of isometries from a metric space to itself forms a group with respect to function composition, called the isometry group.

GOP (Group Of Pictures)

From Wikipedia

In MPEG encoding, a group of pictures, or GOP, specifies the order in which intra-frames and inter frames are arranged.

The GOP is a group of successive pictures within an MPEG-coded video stream. Each MPEG-coded video stream consists of successive GOPs. From the MPEG pictures contained in it the visible frames are generated.

A GOP can contain the following picture types:

I-picture or I-frame (intra coded picture) reference picture, corresponds to a fixed image and is independent of other picture types. Each GOP begins with this type of picture.
P-picture or P-frame (predictive coded picture) contains motion-compensated difference information from the preceding I- or P-frame.
B-picture or B-frame (bidirectionally predictive coded picture) contains difference information from the preceding and following I- or P-frame within a GOP.
D-picture or D-frame (DC direct coded picture) serves the fast advance.

A GOP always begins with an I-frame. Afterwards several P-frames follow, in each case with some frames distance. In the remaining gaps are B-frames. With the next I-frame a new GOP begins.

The GOP structure is often referred by two numbers, for example M=3, N=12. The first one tells the distance between two anchor frames (I or P). The second one tells the distance between two full images (I-frames), it is the GOP length. For the above example, the GOP structure is IBBPBBPBBPBB. Instead of the M parameter one can use the maximal count of B-frames between two consecutive anchor frames.

The more I-frames the MPEG stream has, the more it is editable. However, having more I-frames increases the stream size. In order to save bandwidth and disk space, videos prepared for internet broadcast often have only one I-frame per GOP.

The I-frames contain the full image, they don't require any additional information to reconstruct the image. Therefore any errors in the streams are corrected by the next I-frame (an error in the I-frame propagates until the next I-frame). Errors in the P-frames propagate until the next anchor frame (I or P). B-frames do not propagate errors.

How to Write a Video Player in Less Than 1000 Lines

This is an amazing tutorial of how to use ffmpeg library and SDL to write your OWN video media player~ and it really works! Plus use QT, you may design your own graphic UI as beautiful as you can....

http://www.dranger.com/ffmpeg/

How to find median value of a given two-dimensional array??

Find median is a piece of cake in matlab, but a little bit tricky when comes to C++..

Usually, you need to implement a sort algorithm then find the median as in the middle of the sorted array...

I found the easiest way of doing it is using STL (Standard Template Library).

Here is the method,

for example, you have a two dimensional array

double mydata[mb_width][mb_height];

and you like to find the median value of this array..

To use STL library, you first have to include the header file algorithm.h

#include "algorithm.h"

then define a vector mytheta:

vector mytheta;

pass value from mydata to mytheta;

for (int i=0;i< mb_width;i++){
for (int j=0;j< mb_height;j++){
mytheta.push_back( mydata[i][j]);
}
}

Then you are free to use sort()

sort(mytheta.begin(), mytheta.end());

And here is the median:

median= *(mytheta.begin()+mytheta.size()/2);

Easy, hr?

SDL (Simple Directmedia Layer) Library

Simple DirectMedia Layer is a cross-platform multimedia library designed to provide low level access to audio, keyboard, mouse, joystick, 3D hardware via OpenGL, and 2D video framebuffer. It is used by MPEG playback software, emulators, and many popular games, including the award winning Linux port of "Civilization: Call To Power."

SDL supports Linux, Windows, Windows CE, BeOS, MacOS, Mac OS X, FreeBSD, NetBSD, OpenBSD, BSD/OS, Solaris, IRIX, and QNX. The code contains support for AmigaOS, Dreamcast, Atari, AIX, OSF/Tru64, RISC OS, SymbianOS, and OS/2, but these are not officially supported.

SDL is written in C, but works with C++ natively, and has bindings to several other languages, including Ada, C#, D, Eiffel, Erlang, Euphoria, Guile, Haskell, Java, Lisp, Lua, ML, Objective C, Pascal, Perl, PHP, Pike, Pliant, Python, Ruby, Smalltalk, and Tcl.

What a summary!!

Today I make a summary of all programs i wrote last month, which turned out
to be a long list of files....feels good! :)

*****************************************************************************************

Experimental Summary (Under linux, using gcc, g++, and qt)

*****************************************************************************************

Excutable Program C****** (****** main program)

compiled by using

g++ -o ****** ******.cpp ******.cpp ******.cpp ******.cpp ******.c ******.c -l***** -lav****t

where **********.cpp calls 3 classes,

class ****** has 2 files ( ******.cpp, ******.h)
//class of Camera Motion Estimation, Motion Vector Extraction and drawing on frames,
//saving frames to files.

class ***** has 2 files (*****.cpp, *****.h)
//class of adding captions to frame files, captions include camera motion information

class ** has 2 files (**.cpp, **.h)
//class of estimation minimum of a cost function based on Nelder_Mead algorithm.

-l********: library of *********, from*********
-l********: library of *********, from*********

***.c: define what font is to use, called by ***.
****.c: define ppm file operations, like copy, read, write, called by ****.

*****************************************************************************************

*File_list*

----------------------------------------------------------------

father folder: /**/****/***
includes:

----------------------------------------------------------------
**folders:

CameraPara: folder contains files containing estimated camera parameters.
infile: folder contains files of extracted raw frames from MPEG video.
outfile: folder contains files of output frames, which has motion vectors, as well as camera motion parameter displaying on them.
result bake: bake folder contains past experimental results (Camera para, infile, outfile).
----------------------------------------------------------------
**source files:

*********.cpp: main program to call other programs.
*******.cpp: cpp file for class ****.
***.h: head file for class ****.
***.h: used by class **** for configuration.
******.c: c file define fonts for class ***.
*****.h: head file for *******.c.
********.cpp: test file to test for Camera Parameter Estimation using class
***** and ****. (mannually define camera para);
main1.cpp: test file to test for Camera Parameter Estimation.
****.h: head file called by class ****.
***.bdf: bdf file called by class ******.
*****.cpp: cpp file for class ****.
*****.h: head file for class *****.
******.cpp: same as ******.cpp.
****.h: same as ****.h.
**.c: c file define operations of ppm file.
**.h: head file for **.c
*********.cpp: cpp file for class *********.
*********.h: head file for class **********.
******.h: head file for using stl-vector class
test.cpp: test file used for camera para display on frames.
-----------------------------------------------------------------
**Excutable files:

CMotionEst: main program for camera motion estimation from videos,
#usage: CMotionEst source_mpeg_file.
#compile: g++ -o CMotionEst ******.cpp **.cpp ***.cpp *****.cpp *****.c ****.c -l*****c -l******

nelder: main program for finding a minimum value of function of
multiple variables.
#usage: nelder ( will display a very simple function y=x[0]^2+x[1]^2+x[2]^2+x[3]^2.

ppmaddcaption: main program for testing and display caption on one frame.
#usage: ppmaddcaption

CamEst: main program for testing and display camera motion.
#usage: CamEst
-------------------------------------------------------------------

How to pass a multi-dimensional array as reference to a function??

If you have a multi-dimensional data, e.g. two-dimensional array, say

int mydata[10][20];

and you would like to pass it as reference to another function in C++, say

int myfunc (int data[][])
{
//got mydata[10][20] here..
//do whatever you like to mydata[10][20]
}

How can you do it? well, in C++, it is not easy to handle multi-dimensional array (Weired...why? I DON'T KNOW, but it is the TRUTH!). Well, you may always re-organize the data to a 1D array and simply pass pointer of that array...(which seems easier), but you have to remember the way you put them to 1D, and well, you have to re-organize back to multi-dimensional array after you passed inside the function....I don't prefer that....So, isn't there way that can pass multi-dimensional array as it was to a function?

Well, luckily we have STL, the standard template library. We may define a two-dimensional array by using the vector class in STL, e.g.

> my2Darray(10, vector(20,0));
//this will initialize an 2D array with initial values set to 0;

and then you may pass the 2D array as reference to a function, like this,

void myfunc ( >& anyname)
{
//well, the 2D array "anyname" is already passed inside!
//you may do whatever on it...
//e.g. print out one of its element...
cout<< my2Darray[1][2]<
//firstly define a 2D array, or you may assign other values...

> my2Darray(10, vector(20,0));

//then call myfunc

myfunc(my2Darray);

return 0;

}

You will see the output of print out of my2Darray[1][2];

Bingo!

Dynamic Allocation of Multidimensional Array

I suddenly found that actually the method of dynamic allocation of one-dimensional array is easy, just need a new operation of dynamic allocation of array, but for multi-dimensional array (e.g. matrix, tensor), you CAN'T just use one new operation.

Solutions:

(1) Pre-define a template class which deal with array

#using a standard container to calling
new [].

#include
#include
#include
#include

template <>
struct vector2d
{
vector2d(std::size_t n = 0) : m_Order(n), data( n * n ) {}
void set_Order(std::size_t n) { data.resize(n * (m_Order = n)); }

T *operator[](std::size_t off)
{
return &(data[off * m_Order]);
}
T const *operator[](std::size_t off) const
{
return &(data[off * m_Order]);
}

private:

std::size_t m_Order;
std::vector<> data;
};


int main()
{
using namespace std;

vector2d<> a(10);
int i, j;

for (i = 0; i < j =" 0;" i =" 0;" j =" 0;" style="font-weight: bold;">

(2) Use STL vector class, like std::vector, or std::valarray (recommanded!)

Two / Three / Multi Dimensioned arrays using vector:

A two dimensional array is a vector of vectors. The vector contructor can initialize the length of the array and set the initial value.

Example of a vector of vectors to represent a two dimensional array:

#include 
#include 

using namespace std;

main()
{
// Declare size of two dimensional array and initialize.
vector<> > vI2Matrix(3, vector(2,0));

vI2Matrix[0][0] = 0;
vI2Matrix[0][1] = 1;
vI2Matrix[1][0] = 10;
vI2Matrix[1][1] = 11;
vI2Matrix[2][0] = 20;
vI2Matrix[2][1] = 21;

cout << "Loop by index:" << ii="0;" jj="0;">

A three dimensional vector would be declared as:

#include 
#include 

using namespace std;

main()
{
                           // Vector length of 3 initialized to 0
vector vI1Matrix(3,0);

                           // Vector length of 4 initialized to hold another
                           // vector vI1Matrix which has been initialized to 0
vector<> > vI2Matrix(4, vI1Matrix);

                           // Vector of length 5 containing two dimensional vectors
vector<> > > vI3Matrix(5, vI2Matrix);

...

or declare all in one statement:

#include
#include

using namespace std;

main()
{
vector<> > > vI3Matrix(2, vector<> > (3, vector(4,0)) );

for(int kk=0; kk<4; jj="0;" ii="0;">

Source above (http://www.yolinux.com/TUTORIALS/LinuxTutorialC++STL.html)

(3) define a one-dimensional array of pointers, each points to a one-dimensional array, so that's two dimensional!

(4) convert you problem of multidimensional array to one-dimensional, and remember the way you put them to one-dimensional, then you may do whatever using new operation!

Anyway, multidimensional array is actually one-dimensional array, the ONLY difference is they are put in different format!

Dynamic Allocation of Array with varying sizes

The problems with fixed size arrays

Declaring an array with a fixed size like

int a[100000];

has two typical problems:

• Exceeding maximum. Choosing a real maximum is often impossible because the programmer has no control over the size of the data sets the user is interested in. Erroneous assumptions that a maximum will never be exceeded are the source of many programming bugs. Declaring very large arrays can be extremely wasteful of memory, and if there are many such arrays, may prevent the program from running in some systems.
• No expansion. Using a small size may be more efficient for the typical data set, but prevents the program from running with larger data sets. If array limits are not checked, large data sets will run over the end of an array with disastrous consequences. Fixed size arrays can not expand as needed.

These problems can be avoided by dynamically allocating an array of the right size, or reallocating an array when it needs to expand. Both of these are done by declaring an array as a pointer and using the new operator to allocate memory, and delete to free memory that is no longer needed.

This is exactly what is vector does, but let's see how it's done with an array.

Declare array as a pointer, allocate with new

To create a variable that will point to a dynamically allocated array, declare it as a pointer to the element type. For example,

int* a = NULL;  // pointer to an int, intiallly to nothing.

A dynamically allocated array is declared as a pointer, and must not use the fixed array size declaration. The above declaration creates a pointer, but doesn't yet allocate any memory to it.

Allocate an array with code>new

When the desired size of an array is known, allocate memory for it with the new operator and save the address of that memory in the pointer. Remember: Pointers may be subscripted just as arrays are. The example below reads in a number and allocates that size array.

int* a = NULL;   // Pointer to int, initialize to nothing.

int n;           // Size needed for array

cin >> n;        // Read in the size

a = new int[n];  // Allocate n ints and save ptr in a.

for (int i=0; i   a[i] = 0;    // Initialize all elements to zero.

}

. . .  // Use a as a normal array

delete [] a;  // When done, free memory pointed to by a.

a = NULL;     // Clear a to prevent using invalid memory reference.

Freeing memory with delete

When you are finished with dynamically allocated memory, free it with the delete operator. After memory is freed, it can be reused by later new requests. Memory that your program didn't free will be freed when the program terminates. Never free memory that wasn't dynamically allocated - the results are unpredictable.

delete [] a;  // Free memory allocated for the a array.

a = NULL;     // Be sure the deallocated memory isn't used.

Use [] when deleting arrays

You must specify "[]" when deleting an array, but not for a single value. It isn't possible to delete only part of an array.

Do you have to reset a pointer after delete?

Following the delete in these examples, I reset the pointer to NULL. This isn't strictly necessary, but it's very good practice so that any use of the pointer will produce an error. Attempts to use memory location 0, which is the normal default value of NULL, will be blocked by the way most operating systems allocate memory.

Why doesn't delete reset the pointer? It does in some systems, but the language specification does not require it, so not all systems do it.

Source: http://www.fredosaurus.com/notes-cpp/newdelete/50dynamalloc.html

How to return multiple variables??

In C++, every function will have some return value, however, the number of permitted return value is one.

For example,

int myfunc()
{ int a, b;

a=3; b=4;

return a,b;
}

will cause trouble since the calling function will not know which value is returned...

Then how we return multiple values like described above??

We may define a structure data, for example,

typedef struct Mystruct
{ int a,
int b
} MYSTRUCT;

in the header file,

then in the cpp file,

we define the function as MYSTRUCT type:

MYSTRUCT myfunc()
{
MYSTRUCT mystruct;

mystruct.a=3;
mystruct.b=4;

return mystruct;
}

Here we go!

We may define another MYSTRUCT type data and the returned value would be transfered to that data by the calling function.

Cheers!

What's inline function in c++

When a function is declared inline, the function is expanded at the calling block. The function is not treated as a separate unit like other normal functions.
But a compiler is free to decide, if a function qualifies to be an inline function. If the inline function is found to have larger chunk of code, it will not be treated as an inline function, but as like other normal functions.

Inline functions are treated like macro definitions by the C++ compiler. They are declared with the keyword inline as follows.

//Declaration for C++ Tutorial inline sample:
int add(int x,int y);

//Definition for C++ Tutorial inline sample:
inline int add(int x,int y)
{
return x+y;
}
In fact, the keyword inline is not necessary. If the function is defined with its body directly and the function has a smaller block of code, it will be automatically treated as inline by the compiler.

As implied, inline functions are meant to be used if there is a need to repetitively execute a small block of code, which is smaller. When such functions are treated inline, it might result in a significant performance difference.

WhY Virtual Functions in C++??

C++ virtual function is a member function of a class, whose functionality can be over-ridden in its derived classes. The whole function body can be replaced with a new set of implementation in the derived class. The concept of c++ virtual functions is different from C++ Function overloading.

C++ Virtual Function - Properties:

C++ virtual function is,
A member function of a class
Declared with virtual keyword

Usually has a different functionality in the derived class
A function call is resolved at run-time

The difference between a non-virtual c++ member function and a virtual member function is, the non-virtual member functions are resolved at compile time. This mechanism is called static binding. Where as the c++ virtual member functions are resolved during run-time. This mechanism is known as dynamic binding.

C++ Virtual Function - Reasons:

The most prominent reason why a C++ virtual function will be used is to have a different functionality in the derived class.

For example a Create function in a class Window may have to create a window with white background. But a class called CommandButton derived or inherited from Window, may have to use a gray background and write a caption on the center. The Create function for CommandButton now should have a functionality different from the one at the class called Window.

C++ Virtual function - Example:

This article assumes a base class named Window with a virtual member function named Create. The derived class name will be CommandButton, with our over ridden function Create.

class Window // Base class for C++ virtual function example
{ public:
virtual void Create() // virtual function for C++ virtual function example
{ cout <<"Base class Window"<

class CommandButton : public Window
{ public: void Create()
{ cout<<"Derived class Command Button - Overridden C++ virtual function"<
void main()
{ Window *x, *y;
x = new Window();
x->Create();
y = new CommandButton();
y->Create(); }

The output of the above program will be,
Base class Window
Derived class Command Button

If the function had not been declared virtual, then the base class function would have been called all the times. Because, the function address would have been statically bound during compile time. But now, as the function is declared virtual it is a candidate for run-time linking and the derived class function is being invoked.

C++ Virtual function - Call Mechanism:

Whenever a program has a C++ virtual function declared, a v-table is constructed for the class. The v-table consists of addresses to the virtual functions for classes and pointers to the functions from each of the objects of the derived class. Whenever there is a function call made to the c++ virtual function, the v-table is used to resolve to the function address. This is how the Dynamic binding happens during a virtual function call.

Gauss-Newton algorithm from wikimedia

The Gauss–Newton algorithm is a method used to solve non-linear least squares problems. It can be seen as a modification of Newton's method for finding a minimum of a function. Unlike Newton's method, the Gauss–Newton algorithm can only be used to minimize a sum of squared function values, but it has the advantage that second derivatives, which can be challenging to compute, are not required.
Non-linear least squares problems arise for instance in non-linear regression, where parameters in a model are sought such that the model is in good agreement with available observations.
The method is due to the renowned mathematician Carl Friedrich Gauss.

Conjugate Gradient Methods in Multidimensional non-linear Optimization

Use gradient information:

f(x)~=c-b*x+1/2x*A*x.

Steepest Descent algorithm....not so good.

Solution: Conjugate Gradient Method.

Minimizing the function

f(x)=1/2 *x*A*x-b*x.

The function is minimized when its gradient

gradient(f)=A*x-b is zero.

The conjugate gradient method can be used to solve not only linear algebraic equations by minimizing a quadratic form, but also minimizing a non-linear function.

*Conjugate directions
--"non-interfering" directions with special property that minimizing along one is not "spoiled" by subsequently minimizing along another.

In conjugate gradient method, we want to find new gradient in a direction that is "conjugate" to the old gradient, and to all previous directions traversed.

Powell's Quadratically Convergent Method

Powell first discovered a direction set method that does produce N mutually congugate directions.

Here is how it goes:

Initialize the set of directions u_i to the basis vectors,

u_i=e_i, for i=1,...,N.

Now repeat the following steps until your function stops decreasing:

1. Save your start position as P_0
2. For i=1,..,N. Move P_{i-1} to the minimum along direction u_i and call this point P_i.
3. For i=1,...,N-1, set u_i<-u_{i+1}.
4. Set u_{N}<-P_{N}-P_0
5. Move P_{N} to the minimum along direction u_{N} and call this point P_0.

Done!

Numerical solution to multidimensional non-linear optimization

Basic idea, start at a point P in N-dimensional space, and we proceed from there in some vector direction n, then any function of N variables f(P) can be minimized along the line n by our one-dimensional methods. So a multidimensional minimization methods can be transformed as sequences of line minimizations. Different methods differ only by how, at each stage, they choose the next direction n to try.

The basic routine of line minimization is

linmin: Given as input the vectors P and n, and the function f, find the scalar /lambda that minimizes f(P+\lambda*n). Replace P by P+\lambda*n, replace n by \lambda*n. Done.

G++ and qt in linux

Qt is really nice to develop graphic user interface under linux platform, just the same as in Visual C++..

For me, I prefer gcc+Qt under linux, and Visual C++ under windows.

Well, you can also use Qt under windows, but need to install MSYS and MinGW, which is really lots of extra work! But after that, basically you are able to run a window interface under windows using linux command line, as if you were programming under linux, isn't that cool?

A good book for numerical computation

"Numerical Recipes in C, the art of Scientific Computing", authored by W.H.Press, S.A. Teukolsky, W.T.Vetterling and B.P. Flannery.

I really found this book very useful, since it covers almost all problems in numerical computation, like solution of linear algebraic equations, functions, sorting, minimizing and maximizing of functions, eigensystems..

I read this book since I come to a multidimensional non-linear function minimization problem, and I search on the internet found that there is one very famous algorithm called the "Downhill Simplex Method" by Nelder and Mead.