feval() with fun as a vector

Hi, I am looking how to use feval with the "fun" string input as a vecor

F=[f1,f2,...,fn]

This is in order tu compute the Jacobian

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

function y = test(x)

y=[x(1)^3+x(2)^3 , x(1)];

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

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

function g = jacobian(fun, x0)

delta = 0.001; 
for i = 1 : min(size(x0))
   for j = 1 : max(size(x0))
        x=x0; 
        x(j) = x0(j);
        aux1 = arrayfun ( fun(i), x ); 
        x(j) = x0(j) + delta;
        aux2 = feval ( fun(i), x ); 
        g(i,j)=(aux2-aux1)/(delta);  
   end
end
end
----------------------------------

And I call the function like this

jacobian('test',x)

with x=[1,1]

The problem is in fun(i), that gives a "t" from "test", I would like to take only the "i" element of the vector.

Do you have a clue?

arrayfun and eval do the same

Regards

Answer 1

[Walter Roberson]

Provided that I have understood your requirements...

fbase = str2fun(fun);
ith = @(v,i) v(i);

for i = 1 : min(size(x0))
  fi = @(x) ith(fbase(x), i);
  for j = 1 : max(size(x0))
    [...]
    aux1 = fi(x);
    [...]
    aux2 = fi(x);
    [...]
  end
end

[Julio Cesar]

thank you!Could you find an improvement? (look dowm)

Answer 2

[Julio Cesar]

Thanks to Walter:

function g = jacobian2(fun, x0)

delta = 0.001;

fbase = str2func(fun); ith = @(v,i) v(i);

for i = 1 : max(size(x0))

fi = @(x) ith(fbase(x), i);
for j = 1 : max(size(x0))
  x=x0;
  x(j) = x0(j);
  aux1 = fi(x);
  x(j) = x0(j) + delta;
  aux2 = fi(x);
  g(i,j)=(aux2-aux1)/(delta);
end
end
end

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

function y = test(x) y=[4*x(1)+x(2),6*x(1)*x(2)];

--------------------------- >>x=[2 2]; >>jacobian2('test',x)

ans =

    4.0000    1.0000
   12.0000   12.0000

I followed this steps http://www.monteroespinosa.com/fotos/resueltos/solucion_20406102.gif

[Walter Roberson]

I'm not sure what improvement you were hoping for?

[Walter Roberson]

test = @(x) [4*x(1)+x(2), 6*x(1)*x(2)];

If you were to use that, then you would not need to use str2func.

You could vectorize without much trouble, at least for the test functions you showed. For example,

test = @(x) [4*x(:,1)+x(:,2), 6*x(:,1)*x(:,2)];

Then a single invocation of fi on an N by 2 array of x values would calculate them all. You would, however, need to modify to

ith = @(v,i) v(:,i);

Your statement

x(j) = x0(j);

is seemingly redundant, coming right after

x = x0;

if x has been assigned all of x0, then copying the j'th position of x0 to the j'th position of x is just going to be copying something that is already there.

x1 = x0;

x1(j) = x1(j) + delta;

aux = fi([x0;x1]);

g(i,j) = (aux(2)-aux(1))/delta;

And you can proceed from there to vectorize over the entire j loop. All provided that your test function is vectorizable.