problem with integral2 function in my codes

integral2 function does not work proberly in my case. How can I simplify the expression in int5=integral2(ht,0,(2*L2),0,L1,'RelTol',1e-12,'AbsTol',1e-14,'Method','iterated'); and make the code work? Any help is appreciated.

For details, my codes are:

format long
clear all
close all
clc
syms x
syms y
%syms t

q=0.2; %W/cm2
a=0.25; %cm
L1=0.5; %cm
L2=1.5; %cm
h=0.001; %W/cm2.K
k=0.011; %W/cm.K
T_ambient=23;
T_initial=T_ambient;
alpha=0.0005789;%k/ro*c
delt=0.0220; %delt stands for deltatime<=(delx^2)/2alfa    
Ts=1;

T=0;
for i=1:10
lambda=(i*pi)/(2*L2);
int1=integral(@(x)(q*cos(lambda*x)),L2-a,L2+a);
int2=integral(@(x)((cos(lambda*x)).^2),L2-a,L2+a);
Cn=int1./int2;
A2B1=-Cn./(k*lambda);
p1=(k.*lambda.*cosh(lambda.*L1))+(h.*sinh(lambda.*L1));
p2=(h.*cosh(lambda.*L1))+(k.*lambda.*sinh(lambda.*L1));
B3=p1./p2;
T=T+(A2B1.*cos(lambda.*x)).*(sinh(lambda.*y)-B3.*cosh(lambda.*y));
end
theta=0+T_initial;
for t=1:Ts
    for m=1:5
            mu=(m*pi)/(2*L2);        
        for n=1:5
            beta=-662.09*n+3755.8;
            int3=integral(@(x)(cos(mu.*x).^2),0,2*L2);
            int4=integral(@(y)(cos(beta.*y).^2),0,L1);
            ht=matlabFunction(-T.*cos(mu.*x)*cos(beta.*y));
            int5=integral2(ht,0,(2*L2),0,L1,'RelTol',1e-12,'AbsTol',1e-14,'Method','iterated');
            Fmn=int5/(int3*int4);
            theta=theta+(Fmn.*cos(mu.*x)*cos(beta.*y)*exp(-alpha*(mu.^2+beta.^2).*t));
        end
    end
    figure
    ezcontourf(theta+T,[0 3],[0 0.5]);
    xlabel('x [cm]')
    ylabel('y [cm]')
    title('Temperature Distribution [C]')
end

Answer 1

[Star Strider]

Your code works (in R2016b). The integral2 function is encountering a discontinuity, and throwing the appropriate warnings.

Plotting your ‘ht’ function will provide you with insight into the problem:

ht = @(x,y)cos(x.*pi.*(1.0./3.0)).*cos(y.*2.43162e3).*(cos(x.*pi.*(5.0./3.0)).*(cosh(y.*pi.*(5.0./3.0)).*1.010332969266712-sinh(y.*pi.*(5.0./3.0))).*6.32895191273502e-16-cos(x.*pi.*(8.0./3.0)).*(cosh(y.*pi.*(8.0./3.0)).*1.000450157978233-sinh(y.*pi.*(8.0./3.0))).*2.262610130858262-cos(x.*pi.*(4.0./3.0)).*(cosh(y.*pi.*(4.0./3.0)).*1.029468631679565-sinh(y.*pi.*(4.0./3.0))).*5.079090139740538+cos(x.*pi).*(cosh(y.*pi).*1.085034525163738-sinh(y.*pi)).*5.52568847333794e-15+cos(x.*pi.*(1.0e1./3.0)).*(cosh(y.*pi.*(1.0e1./3.0)).*1.000055664759626-sinh(y.*pi.*(1.0e1./3.0))).*7.94622746325111e-1+cos(x.*pi.*(7.0./3.0)).*(cosh(y.*pi.*(7.0./3.0)).*1.001279352938934-sinh(y.*pi.*(7.0./3.0))).*2.290646007531096e-15+cos(x.*pi.*(2.0./3.0)).*(cosh(y.*pi.*(2.0./3.0)).*1.254534770465853-sinh(y.*pi.*(2.0./3.0))).*9.074926360779978-cos(x.*pi.*(1.0./3.0)).*(cosh(y.*pi.*(1.0./3.0)).*1.836310670638647-sinh(y.*pi.*(1.0./3.0))).*1.395001209824949e-13-cos(x.*pi.*3.0).*(cosh(y.*pi.*3.0).*1.000158327683549-sinh(y.*pi.*3.0)).*5.797477778816949e-16+cos(x.*pi.*2.0).*(cosh(y.*pi.*2.0).*1.003634943955502-sinh(y.*pi.*2.0)).*3.684406677903193);

x = linspace(0, 3.0, 50);
y = linspace(0, 0.5, 50);
[X,Y] = meshgrid(x,y);

figure(2)
meshc(X, Y, ht(X, Y))
grid on