Prev: ReplaceAll reloaded
Next: Clip[DirectedInfinity[1]] remains unevaluated?
From: Kamil on 11 Feb 2010 06:12
Dear all,

I am kind of new to this forum and mathematica but I have found the
post here quite interesting and helpful. I am having problem
translating this algorithm to mathematica code and I hope I can be
helped. The algorith is as given below. Thanks to you all.

Suppose that for a given matrix A ( n x n ) that is a function of say
x the non-trivial solution is desired.
The following algorithm is used to compute the value of x that makes
the determinant of A zero.

INPUT
(i) The matrix A and find matrix B, where B = - dA/dx (i.e derivative
of A wrt x)
(ii) The initial estimate of x
(iii) A tolerance for convergence say delta =EF=80 (a small positive number)
is selected.
ITERATION
(a) Choose the initial guess x=EF=80 (0) and start the iteration.
(b) Compute the eigenvalues of the matrix (A - B) based upon the
initial guess x=EF=80 (0) .
(c) Evaluate the minimum eigenvalue of found in step b and assign this
value to epsilon
(d) Compute the new estimate x=EF=80 (1)=x(0) + epsilon
(e) Compute the matrices A(x) and B(x) by substituting x= x(1)

(f) Repeat steps (b)-(e) for kth iteration until the condition
Abs(epsilon.... in step c) < delta in step (iii) =EF=80 is satisfied.
(g) Stop the iteration.
(h) Store the value of x(k)
(i) Repeat steps (a)-(h) to evaluate another x(k) =EF=80 for different
starting value x=EF=80 (0) .

From: Daniel Lichtblau on 12 Feb 2010 04:39
Kamil wrote:
> Dear all,
>
> I am kind of new to this forum and mathematica but I have found the
> post here quite interesting and helpful. I am having problem
> translating this algorithm to mathematica code and I hope I can be
> helped. The algorith is as given below. Thanks to you all.
>
> Suppose that for a given matrix A ( n x n ) that is a function of say
> x the non-trivial solution is desired.
> The following algorithm is used to compute the value of x that makes
> the determinant of A zero.
>
> INPUT
> (i) The matrix A and find matrix B, where B = - dA/dx (i.e derivative
> of A wrt x)
> (ii) The initial estimate of x
> (iii) A tolerance for convergence say delta =EF=80 (a small positive number)
> is selected.
> ITERATION
> (a) Choose the initial guess x=EF=80 (0) and start the iteration.
> (b) Compute the eigenvalues of the matrix (A - B) based upon the
> initial guess x=EF=80 (0) .
> (c) Evaluate the minimum eigenvalue of found in step b and assign this
> value to epsilon
> (d) Compute the new estimate x=EF=80 (1)=x(0) + epsilon
> (e) Compute the matrices A(x) and B(x) by substituting x= x(1)
>
> (f) Repeat steps (b)-(e) for kth iteration until the condition
> Abs(epsilon.... in step c) < delta in step (iii) =EF=80 is satisfied.
> (g) Stop the iteration.
> (h) Store the value of x(k)
> (i) Repeat steps (a)-(h) to evaluate another x(k) =EF=80 for different
> starting value x=EF=80 (0) .

If you are not tied to the particular method you indicate above, you
could just use FindRoot on a suitably defined numeric function.

numericalDet[mat:{{_?NumericQ..}..} /;
MatrixQ[mat] && Length[mat]==Length[mat[[1]]]] := Det[mat]

Here is an example.

matrix[x_] := {
{Cos[x], Sin[x], x},
{x-1, x+3/7, Cos[x]^2},
{-x^2+1/3, Sin[x]*Cos[x], 1}}

In[29]:= detroot = FindRoot[numericalDet[matrix[x]], {x,-2}]
Out[29]= {x -> -1.47929}

Check the resulting determinant:

In[30]:= numericalDet[matrix[x]/.detroot] // InputForm
Out[30]//InputForm= -1.628393168080853*^-15

Daniel Lichtblau
Wolfram Research

 | 
Pages: 1
Prev: ReplaceAll reloaded
Next: Clip[DirectedInfinity[1]] remains unevaluated?