How can I integrate this function to give a symbolic expression?

1 view (last 30 days)
I am trying to compute the above double integral to give a function C which depends on x. I attempted to use the following code to compute the double integral: (c,a,b,Da,Beta and sigma are known constants and C_s(x,t) is known to be exp(x*t) in this example. I am trying to plot the function which results with C on the y-axis and x on the x-axis. What is wrong with my integral?
c=0;
Da = .4
sigma = .1
a=1
b=1.3
L=Inf
t=1
Greens= @(x,zeta,t) (1/2*sqrt(pi*a*t))*exp(b(zeta-x)/2*a + c-b^2/4*a)*(exp(-(x-zeta)^2/4*a*t)+exp(-(x+zeta)^2/4*a*t))
fun = @(x,t,zeta,tau) Da.*sigma.*exp(x.*t).*(1/2.*sqrt(pi.*a.*t)).*exp(b(zeta-x)/2.*a + c-b^2/4.*a).*(exp(-(x-zeta)^2/4*a*t)+exp(-(x+zeta)^2/4*a*t))
upper=@(tau) t
q=integral2(fun,0,L,0,upper)
I am getting the following error: "Error using arrayfun All of the input arguments must be of the same size and shape. Previous inputs had size 45 in dimension 2. Input #4 has size 1"

Accepted Answer

Star Strider
Star Strider on 11 Jun 2018
It will not give you a symbolic expression because you are not using the Symbolic Math Toolbox to integrate it. (That would likely not work anyway.) You can have it produce an anonymous function:
c=0;
Da = .4;
sigma = 0.1;
a=1;
b=1.3;
L=Inf;
Greens= @(x,zeta,t) (1/2*sqrt(pi*a*t)).*exp(b*(zeta-x)/(2*a) + (c-b^2/(4*a))*t).*(exp(-(x-zeta).^2/(4*a*t))+exp(-(x+zeta).^2./(4*a*t)));
fun = @(x,t,zeta,tau) Da.*sigma.*exp(x.*t).*Greens(x,zeta,t-tau);
C = @(x,t) integral2(@(zeta,tau)fun(x,t,zeta,tau),0,t,0,L);
X = linspace(0,0.05,10);
T = linspace(0,0.08,20);
for k1 = 1:numel(X)
for k2 = 1:numel(T)
Z(k1,k2) = C(X(k1),T(k2));
end
end
figure
surf(T, X, Z)
There were several errors in your code that I corrected. This runs, although it does not give any useful output (all are NaN). I leave you to solve that problem, and make any corrections that may be necessary, since I have no idea what you are doing or what the correct values for ‘x’ and ‘t’ are.

More Answers (0)

Products


Release

R2017a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!