Integration of these kind of functions

4 次查看(过去 30 天)
I want to integrate
is a 2 by 2 matrix (given below), is identity matrix of order 2, and Eand η and positive real numbers. I want to calculate this integration (for each component of the matrix) with .
Note that the integrand does not diverges as long as η is non-zero, however when η is significantly small, the integrand can be very large for some points in plane. As an example, I have plot the imaginary part of (1,2) componet of the integrand matrix for and below. (the other components also follow the same behavior):
It can be seen that when η is small the integrand is fairly smooth and can be integrated easily even with integrate2 function. But when η is small, the integrand has very sharp peaks due to which integral2 fails. To calculate the correct integration with small η by trapz, we need a lot of grid points which makes everything very slow.
Is there a way to calculate this integration with ?
I have analytically calculated the condition for which the integrand shows sharp peaks, that condition is:
where are components of matrix H. We could also find the exact points on which the condition is satisfied.
where x are the roots of equation and y is all values between -1 and +1. And
So, for each y there will be 4 roots, x, which will give eight pionts (four and four . From these eight pionts, we will take only the real ; imaginary must be neglected. Also for each y there will be two points (one and one .
Is there any way to tell MATLAB to take a lot of gird pionts near points and then calculate the integration? Or is there any way to calculate this type of integration?
Code:
clear; clc;
% parameters
E = 4.4;
eta = 0.1;5e-4;
% kx and ky limits and points
dkx = 0.006;
dky = dkx;
xmin = -2*pi/(3*sqrt(3));
xmax = 4*pi/(3*sqrt(3));
ymin = -2*pi/3;
ymax = 2*pi/3;
kxs = xmin:dkx:xmax;
kys = ymin:dky:ymax;
NBZx = length(kxs);
NBZy = length(kys);
% H matrix:
J = 1;
Dz = 0.5;
S = 1;
s0 = eye(2);
sx = [0, 1; 1, 0];
sy = [0, -1i; 1i, 0];
sz = [1, 0; 0, -1];
h0 = 3 * J * S;
hx = @(kx, ky) -J * S * (cos(ky / 2 - (3^(1/2) * kx) / 2) + cos(ky / 2 + (3^(1/2) * kx) / 2) + cos(ky));
hy = @(kx, ky) -J * S * (sin(ky / 2 - (3^(1/2) * kx) / 2) + sin(ky / 2 + (3^(1/2) * kx) / 2) - sin(ky));
hz = @(kx, ky) -2 * Dz * S * (sin(3^(1/2) * kx) + sin((3 * ky) / 2 - (3^(1/2) * kx) / 2) - sin((3 * ky) / 2 + (3^(1/2) * kx) / 2));
H = @(kx, ky) s0 * h0 + sx * hx(kx, ky) + sy * hy(kx, ky) + sz * hz(kx, ky);
%integrand:
G00 = @(kx, ky) inv(E*eye(2) - H(kx,ky) + 1i*eye(2)*eta); %integrand
  2 个评论
David Goodmanson
David Goodmanson 2024-1-17
Hi Luqman,
since G00 appears to be a 2x2 matrix, what quantity are you plotting?
Luqman Saleem
Luqman Saleem 2024-1-17
@David Goodmanson oops, I forget to mention that. It is the imaginary part of (1,2) componet of G00.

请先登录,再进行评论。

回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Numerical Integration and Differentiation 的更多信息

产品


版本

R2023b

Community Treasure Hunt

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

Start Hunting!

Translated by