quadgk problems

The complete problem is very complex, and as in my last question regarding the same project I skipped the introduction, and I will do so again here.

The specific problem is (this time) an Array (2 cells) of function handles, and I am getting the same error every time (regardless of the changes I make in it (* -> .*, A'*B -> (A*B)' and so on...). So here is the code:

z = sym ('z','real');
zd = sym ('zd','real');

%for demonstration the values of M,w,mi,eps... are simplified
M=4;mi=1;eps0=1;w=100;a=1;b=1;h=10;
deltaz=h/M;

%functions g, R & k
k = w * sqrt(mi * eps0);
R = @(z,zd) sqrt((z-zd).*(z-zd) + a.*a);
g = @(z,zd) (exp(-1i.*k.*R(z,zd))./R(z,zd));

for r = 1:1:M 
    for s = 1:1:M

          %N for index s
          Ni1 = @(z) ((s + 1)*deltaz - z)./deltaz;
          Ni2 = @(z) (z - s*deltaz)./deltaz;
          Ni = {Ni1 Ni2}';
          %N for index j
          Nj1 = @(z) ((r + 1)*deltaz - z)./deltaz;
          Nj2 = @(z) (z - r*deltaz)./deltaz;
          Nj = {Nj1 Nj2}';
          %N'
          Npi1 = @(zd) ((s + 1)*deltaz - zd)./deltaz;
          Npi2 = @(zd) (zd - s*deltaz)./deltaz;
          Npi = {Npi1 Npi2}';

          %D = d(N) / d(z)
          D1 = -1./deltaz;
          D2 = 1./deltaz;
          D = cell2mat({D1 D2})';
          %D'
          Dp1 = D1;
          Dp2 = D2;
          Dp = cell2mat({Dp1 Dp2})';
          %FIRST PROBLEM
          funkcija1 = @(z,zd)(Dp*g(z,zd))'; 
          funkcija2 = @(z)D*quadgk(@(zd)funkcija1(z,zd),s*deltaz,(s+1)*deltaz);
          IntPrvi = quadgk(funkcija2(z),r*deltaz, (r+1)*deltaz);

          %SECOND PROBLEM
          funk3 = @(z,zd) (Npi(zd)*g(z,zd))';
          funk4 = @(z)Nj(z)*int(funk3,zd,s*deltaz, (s+1)*deltaz);
          IntDrugi = k.*k.*quadgk(@(z)funk4(z), r*deltaz, (r+1)*deltaz);

          %here comes a blok of code that uses IntPrvi & IntDrugi
      end
  end

for the first integral the problem is matrix sizes and return value type; TRACE:

_??? Error using ==> quadgk>finalInputChecks at 476 Supported classes are 'double' and 'single'.

Error in ==> quadgk>evalFun at 358 finalInputChecks(x,fx);

Error in ==> quadgk>f1 at 375 [y,too_close] = evalFun(tt);

Error in ==> quadgk>vadapt at 269 [fx,too_close] = f(x);

Error in ==> quadgk at 208 [q,errbnd] = vadapt(@f1,interval);

Error in ==> struja>@(z)D*quadgk(@(zd)funkcija1(z,zd),s*deltaz,(s+1)*deltaz) at 94 funkcija2 = @(z)D*quadgk(@(zd)funkcija1(z,zd),s*deltaz,(s+1)*deltaz); %D*int(funk1)

Error in ==> struja at 95 IntPrvi = quadgk(funkcija2(z),r*deltaz, (r+1)*deltaz); _ ??? Error using ==> quadgk>finalInputChecks at 476 Supported classes are 'double' and 'single'.

Error in ==> quadgk>evalFun at 358 finalInputChecks(x,fx);

Error in ==> quadgk>f1 at 375 [y,too_close] = evalFun(tt);

Error in ==> quadgk>vadapt at 269 [fx,too_close] = f(x);

Error in ==> quadgk at 208 [q,errbnd] = vadapt(@f1,interval);

Error in ==> struja>@(z)D*quadgk(@(zd)funkcija1(z,zd),s*deltaz,(s+1)*deltaz) at 94 funkcija2 = @(z)D*quadgk(@(zd)funkcija1(z,zd),s*deltaz,(s+1)*deltaz); %D*int(funk1)

Error in ==> struja at 95 IntPrvi = quadgk(funkcija2(z),r*deltaz, (r+1)*deltaz);_

for the second integral the problem is how to call the Function handle Array Ni/Nj with the variable z, and multiply it (matrix multipl.) with g(z,zd) or with the result of quadprog for funk4. TRACE: _ ??? Subscript indices must either be real positive integers or logicals.

Error in ==> struja>@(z)Nj(z)*int(funk3,zd,s*deltaz,(s+1)*deltaz)

Error in ==> struja>@(z)funk4(z) at 100 IntDrugi = k.*k.*quadgk(@(z)funk4(z), r*deltaz, (r+1)*deltaz);

Error in ==> quadgk>evalFun at 357 fx = FUN(x);

Error in ==> quadgk>f1 at 375 [y,too_close] = evalFun(tt);

Error in ==> quadgk>vadapt at 269 [fx,too_close] = f(x);

Error in ==> quadgk at 208 [q,errbnd] = vadapt(@f1,interval);

Error in ==> struja at 100 IntDrugi = k.*k.*quadgk(@(z)funk4(z), r*deltaz, (r+1)*deltaz); _

I would be VERY gratefull for any insight why these errors are popping up, or any suggestion for a solution

Answer 1

[Tomislav Broni?]

i solved the problem by dividing it up in 4 integrals, since the N vector isnt big (2×1), so I calculated all 4 elements of the 2×2 matrix and then put them back together into a matrix

Answer 2

[Walter Roberson]

For your second problem, your line

funk4 = @(z)Nj(z)*int(funk3,zd,s*deltaz, (s+1)*deltaz);

is one of the few places that does not bind both z and zd as the dummy argument names for the anonymous function. The zd referred to will thus be the sym() that you created with that name.

The int() of a function is a symbolic operation and will always return a symbolic value, so funk4 will return a symbolic value, which quadgk() will not be able to deal with. If you are sure that the integral will be fully resolved to a numeric value, you could wrap the int() in double() to get the MATLAB floating point equivalent of the symbolic number. If you are expecting the int() to return a symbolic function, you may need to use matlabFunction() to convert that symbolic expression in to an anonymous function that quadgk() could integrate over.

I have not worked through your first problem, but your line

IntPrvi = quadgk(funkcija2(z),r*deltaz, (r+1)*deltaz);

is going to be working with symbolic expressions it appears to me.

Answer 3

[Tomislav Broni?]

I already tried all sorts of combinations with the syntax...

dblquad cant work, becouse I need to integrate once, then multiply it with N (or D), and then integrate again.

I have to admit that I'm currently trying out quadgk syntax and i cant get it to work with matrices or vectors of inputs...

I'll try the conversion to a function from sym, but that is just one of the ways I was trying to get it to work.

for both IntPrvi and IntDrugi I have to calculate integral(2×1 * integral(1×2 * g(z,zd))) where in the second the matrices N contain function handles in cells, which is also a great problem to get working...

Thanks for the help, i really could use some more.