计算科学 - 在 MATLAB 中解决冷凝物密度问题 - 吾爱随笔录

我想解决 $n_{0}$ 对于固定值 $n$ ，让我们说 $n=1$

n = n_{0} + \frac{1}{2} \int_{- 1 / 2}^{1 / 2} d q (\frac{e_{q} + U n_{0}}{ℏ ω} - 1)

$n= n_{0}+ \dfrac{1}{2}\int_{-1/2}^{1/2}dq\left(\dfrac{e_{q}+Un_{0}}{\hbar\omega}-1\right)$ 在哪里

e_{q} = 2 [1 - c o s (2 π q)]

$e_{q}=2[1-cos(2\pi q)]$ 和

ℏ ω = (e_{q}^{2} + 2 U n_{0} e_{q})^{1 / 2}

$\hbar\omega=(e_{q}^{2}+2Un_{0}e_{q})^{1/2}$

然后我想在 $U$ 和 $n_{0}$ . 以下是MATLAB中的代码。

 eq = @(q) 2*(1-cos(2*pi*q));
 hq = @(q,U,n0) ((eq(q))^2+2*U*n0*(eq(q))).^(1/2);
 y = @(q,U,n0) (((eq(q))+(U*n0))/hq(q,U,n0))-1;
 a = -0.5;
 b = 0.5;
 v = @(U,n0) quad(@(q) y(q,U,n0),a,b);
 cv =@(U,n0) n0+(0.5*v)-1;
 n0 = 0.1;
 options = optimset('Display', 'iter'); 
 n = @(U) fsolve(vv,n0,0.1);
 U = 0:0.1:20;
 plot(U,n(U))

编辑 2 代码更正编辑 1。

编辑 3

但我仍然面临以下错误：

错误使用 ^ 输入必须是标量和方阵。要计算元素 POWER，请改用 POWER (.^)。

@(q,n0)((eq(q))^2+2*U n0 (eq(q))) 中的错误。^0.5

@(q,n0)(((eq(q))+(U*n0))/hq(q,n0))-1 中的错误

@(q)y(q,n0) 中的错误

quad 中的错误（第 72 行）y = f(x, varargin{:});

@(n0)quad(@(q)y(q,n0),a,b) 中的错误

@(n0)n0+(0.5*v(n0))-1 中的错误

simulv1 中的错误（第 12 行） n = fsolve(cv(n0),n0,0.1,options);