help-octave
[Top][All Lists]

## matrix problems

 From: jonnysukes Subject: matrix problems Date: Tue, 5 Feb 2008 08:48:36 -0800 (PST)

```I'm trying to run a matlab script for a class I'm taking now, but I can't get
it to work. I already went through the script and changed all the syntax to
follow octave notations. The problem I'm having is trying to run the
operation:

C =D\f

I've gone through all the checks in the code and I found that the matrix D
has entries of +-Inf and so it can't run the operation. I'm almost positive
that's the problem because thats where the code stops running. But the
entries in D aren't actually Inf, they're just bigger than double precision
supports. Is there anyway around that problem, or an operation similar to \
that works on infinity? Here's the code as of now, I tried fixing the
problem by just replacing Inf with 2*2**1022, but it doesnt seem to be
working well. I also tried using the symbolic package, but it doesn't seem
to let you index matrices. Thanks for the help!

function nmse_error = refine_prony_hildebrand_error(ri, z); %(ri, t, f)

t = z{1};

f = z{2};

N = length(t);
n = length(ri)/2;

delta = t(2) - t(1);
x = t/delta;

a = ri(1:n) + i*ri((n+1):end);

% compute mu

mu = exp(a);

% construct matrix equation for series coefficients

D = ones(N,n);
bigconstant = 2*2**1022;
for r = 1:N
for s = 1:n
b = mu(s)^(r-1);
if (abs(real(b)) == Inf && abs(imag(b)) == Inf)
D(r,s) = sign(real(b))*bigconstant + sign(imag(b))*i*bigconstant;
elseif(abs(real(b)) == Inf)
D(r,s) = sign(real(b))*bigconstant + i*imag(b);
elseif(abs(imag(b)) == Inf)
D(r,s) = real(b) + sign(imag(b))*i*bigconstant;
else
D(r,s) = mu(s)^(r-1);
end
end
end

% solve system

C = D\f;

% compute Prony fit

f_prony = zeros(size(t));
for k = 1:length(a)
f_prony  = f_prony + exp(a(k)*x)*C(k);
end

% compute error

nmse_error = this_nmse(f, f_prony)

plot(t, real(f), t, real(f_prony))
drawnow

%
% computes the Normalized Mean Square Error between two vectors
%
% e = nmse(f,g)
%

function e = this_nmse(f,g)
e = sum(sum((abs(f-g)).^2)/sum((abs(f)).^2));
--
View this message in context:
http://www.nabble.com/matrix-problems-tp15292917p15292917.html
Sent from the Octave - General mailing list archive at Nabble.com.

```