用matlab求欧拉常数代码，谢谢各位

如何用matlab求欧拉常数？

1、首先我们根据欧拉常数的定义，写出其表达式，如下图所示。

2、从表达式我们看到，求和部分可以用symsum函数来求解

3、然后再用limit函数，求其n一﹥∞的极限

4、完整的代码如下

>>syms k n

>>S = symsum(1/k,k,1,n) - log(n)

>>vpa(limit(S,n,Inf),20)

5、也可以直接用下列命令来求解

>>-psi(1)

6、执行结果

% 这一类计算可以用matlab的符号计算实现

% 首先是定义符号

syms x y z r a b c;

%%问题一：计算极限 采用limit函数

fprintf("---计算极限---\n")

fun1 = ((1+tan(x)) / (1+sin(x)))^(1/x^3);

fprintf("x ->0, fun1->%s\n", limit(fun1));

fun2 = x y / sqrt(x^2 + y^2);

fprintf("(x, y)->0, fun2->%s\n", limit(limit(fun2, x, 0), y, 0));

%% 问题二：计算导数 采用diff函数

fprintf("---计算导数---")

f = exp(-2x) cos(3 sqrt(x));

fprintf("\nf的导数:\n %s\n", diff(f, x));

%% 问题三: 计算积分 采用int函数

fprintf("---积分---");

fun1 = sin(xy + z);

fprintf("\nfun1对z积分: %s\n", int(fun1, z));

fun2 = 1 / (3 + 2x + x^2);

fprintf("fun2对x从0到1积分:%s\n", int(fun2, x, [0, 1]));

fun3 = 1 / (3 + 2x + x^2);

fprintf("fun3对x从-inf到inf积分: %s\n", int(fun3, x, [-inf, inf]));

%% 计算Jacobi矩阵

% 需要知道(x, y, z)是行向量还是列向量

% 假设是行向量

fprintf("--jacobi--");

vector = [rsin(a)cos(b), rsin(a)sin(b), rcos(a)];

fprintf("\nJacobi矩阵\n")

jacobi = [diff(vector, r);

diff(vector, a);

diff(vector, b)];

disp(jacobi)

% 如果是列向量, 只需要对jacobi做一个转置 jacobi=jacobi';

%% 问题五: 证明等式 需要用simplify化简函数

fprintf("--等式证明--")

left = int(x^3cos(ax)^2, x);

right = x^4/8 + (x^3/4/a - 3x/8/a^3)sin(2ax) + (3x^2/8/a^2 - 3/16/a^4)cos(2ax)+c;

% 只需证明两者的一阶导数相同即可 (因为积分允许差一个常数项)

fprintf("\n左右两边一阶导数的差为:%s\n", simplify(diff(left, x) - diff(right, x)));

-----------------------以下是在我matlab中的测试结果---------------------

% 950944444 (R2018b)

---计算极限---

x ->0, fun1->exp(1/2)

(x, y)->0, fun2->0

---计算导数---

f的导数:

- 2exp(-2x)cos(3x^(1/2)) - (3exp(-2x)sin(3x^(1/2)))/(2x^(1/2))

---积分---

fun1对z积分: -cos(z + xy)

fun2对x从0到1积分:(2^(1/2)atan(2^(1/2)/4))/2

fun3对x从-inf到inf积分: (2^(1/2)pi)/2

--jacobi--

Jacobi矩阵

[ cos(b)sin(a), sin(a)sin(b), cos(a)]

[ rcos(a)cos(b), rcos(a)sin(b), -rsin(a)]

[ -rsin(a)sin(b), rcos(b)sin(a), 0]

--等式证明--

左右两边一阶导数的差为:0

如果您学过极限应该对这个函数不陌生。F是一个用x表示的符号表达式，LIMIT（F，x,a）表示当x趋于a时F的极限。

如果F表达式中只有一个变量，那么就没有必要指出x了，直接用LIMIT（F，a）即可。

以此类推。LIMIT（F） 此时默认a = 0相当于LIMIT（F，0）

极限有左右极限之分，因此 可以用‘left’，‘right’分别表示左右极限。

limit求极限，例如1/(x^2)x趋于 负无穷时的值为0

syms x

>> limit(1/(x^2),x,-inf)

ans =

0