Problems with quad

I'm trying to get a numerical value from a function using quad. In a separate .m-file I got a function looking like this:

function Wi = Wifun(y,T,TM,TMO,G2,A,b,x,Epsi,ny)

syms G2;syms T;syms TM;syms TMO;syms ny;syms A;syms Epsi;syms b;syms x;syms y;

G = G2.*(1+((T-300)./TM).*TMO);
Ei = -(1/2.*pi).*((1+ny)./(1-ny)).*G.*A.*Epsi;
Wi = -(2/3).*b.*x.*Ei.*(1./(x.^2+y.^2));

This function will produce this formula:

-(pi*A*Epsi*G2*b*x*((TMO*(T - 300))/TM + 1)*(ny + 1))/(3*(x^2 + y^2)*(ny - 1))

If I write like this:

syms y
T = 1273;
TM = 2163;
TMO = -0.5;
G2 = 126000;
A = 1.2E-29;
b = 0.000000000258;
x = 2/3*b;
Epsi = 0.00487903650119246;
ny = 0.3;

fun = @(y)-(pi.*A.*Epsi.*G2.*b.*x.*((TMO.*(T - 300))./TM + 1).*(ny + 1))./(3.*(x.^2 + y.^2).*(ny - 1));
z2 = quad(fun,-2,2)

It works and z2 gets the value of 4.0641e-027. But if instead of the 2 last lines try to write like this:

z = quad(Wifun,-2,2)

There will be this error:

??? Error using ==> fcnchk at 108
If FUN is a MATLAB object, it must have an feval method.

Error in ==> quad at 66
f = fcnchk(funfcn);

Error in ==> test at 19
z = quad(Wifun,-2,2)

Someone told me that you should write like this:

Wifun2 = str2func('Wifun');
z = quad(Wifun2,-2,2)

But then the error looks like this instead:

??? Undefined function or method 'isfinite' for input arguments of type 'sym'.

Error in ==> quad at 81
if ~isfinite(y(1))

Error in ==> test at 19
z = quad(Wifun2,-2,2)

Is there someone here that knows how I can do this?

[Andrew Newell]

Could you please format the code? See http://www.mathworks.com/matlabcentral/answers/help/markup.

[Jon]

Sorry, didn't know you could do that. Thanks for telling me. Seems like someone already did it for me, thanks :).

Answer 1

[Walter Roberson]

Use matlabFunction() to convert the symbolic expression in to an executable routine.

That executable routine will want to take all those names as parameters but the routine you pass to quad() should only accept one parameter. To take care of this, follow the documentation about "Passing additional parameters" to quad and so on. Which, in short, would look like

wifun_handle = matlabFunction(wifun, {'y' 
T' 'TM' 'TMO' 'G2' 'A' 'b' 'x' 'Epsi' 'ny'});
fun = @(y) wifun_handle(y,T,TM,TMO,G2,A,b,x,Epsi,ny);

and then pass fun to quad.

By the way, you can change your line

function Wi = Wifun(y,T,TM,TMO,G2,A,b,x,Epsi,ny)

to

function Wi = Wifun

as you never actually pass any parameters to Wifun itself.

Answer 2

[Jon]

When trying matlabfunction() I get another problem.

Wi2 = matlabFunction(Wifun)

Produces the following:

Wi2 = 

    @(A,Epsi,G2,T,TM,TMO,b,ny,x,y)(A.*Epsi.*G2.*pi.*b.*x.*((TMO.*(T-3.0e2))./TM+1.0).*(ny+1.0).*(-1.0./3.0))./((x.^2+y.^2).*(ny-1.0))

Then:

fun = @(y) Wi2(y,T,TM,TMO,G2,A,b,x,Epsi,ny);
z = quad(fun,-2,2)

Produces the value z = -1.3553e-020. If I instead write the function manually and solve it with quad like this:

fun3 = @(y) (A.*Epsi.*G2.*pi.*b.*x.*((TMO.*(T-3.0e2))./TM+1.0).*(ny+1.0).*(-1.0./3.0))./((x.^2+y.^2).*(ny-1.0))
z3 = quad(fun3,-2,2)

It produces the value z3 = 4.0641e-027 which is clearly different from the one above.

[Walter Roberson]

You are not passing the parameters to Wi2 in the same order. In particular notice that Wi2 expects Epsi as the second parameter, but your fun passes it as the second-last parameter.

The form of matlabFunction() that I showed is not the one that you used. I showed passing in the list of variable names; you did not. When you do not pass in the list of variable names, then matlabFunction expects the parameters in alphabetical order.

Answer 3

[Jon]

When I do as you wrote:

Wi2 = matlabFunction(Wifun, {'y' 'T' 'TM' 'TMO' 'G2' 'A' 'b' 'x' 'Epsi' 'ny'})
fun = @(y) Wi2(y,T,TM,TMO,G2,A,b,x,Epsi,ny)
z = quad(fun,-2,2)

I get the error:

??? Error using ==> deal at 38
The number of outputs should match the number of inputs.

Error in ==>
sym.matlabFunction>makeFhandle/@(A,Epsi,G2,T,TM,TMO,b,ny,x,y)deal((A.*Epsi.*G2.*pi.*b.*x.*((TMO.*(T-3.0e2))./TM+1.0).*(ny+1.0).*(-1.0./3.0))./((x.^2+y.^2).*(ny-1.0)),[y,T,TM,TMO,G2,A,b,x,Epsi,ny])

Error in ==> @(y)Wi2(y,T,TM,TMO,G2,A,b,x,Epsi,ny)

Error in ==> quad at 77
y = f(x, varargin{:});

Error in ==> test at 20
z = quad(fun,-2,2)

This is why I tried another way.

[Walter Roberson]

Looks like you might have to insert a keyword,

matlabFunction(Wifun, 'Vars', {'y' 'T' <etc>})

(Sorry; I don't have the toolbox myself so I work from the documentation and miss parameters sometimes.)

[Jon]

No problems. Now it is working :). Thank you.