一种MATLAB中解复杂方程（高次、指数、无解析解）的方法，可以在实现论文中公式时使用，solve函数。

      前几天我在实现一篇论文时，对于一个公式其他参数都已知的情况下，要得到剩下得那个未知的变量，由于方程的形式很复杂，用常规方法很难处理，故在实现时使用了MATLAB中solve函数，现在把方法呈现在这里，在复现论文公式时是一种比较简洁的思路。

       首先，我找了一个简单的例子：

σ02-Eγ22l224σ022=σ01-Eγ12l224σ012-αE(t2-t1)

        此公式中，除σ02外其余的参数皆可知，要求输入其余变量，输出为σ02。

        可见，只看σ02情况下，该方程为一三次方程，已经很难用笔算出。那么不妨将其余变量视为已知量，在MATLAB中使用root函数，找到σ02对应其余变量的表达式。

代码：

syms x m n;

f1 = [x^3-m*x^2-(1/24)*n==0];

solve(f1,x)

       由于参数过多，syms对于太多变量时会报错，可以取巧一点用将多个已知变量放在一起用一个变量表示，等号右边用m表示，等号左边第二项分子用n表示，σ02用x表示。

运行结果：

root(z^3 - m*z^2 - n/24, z, 1)
 root(z^3 - m*z^2 - n/24, z, 2)
 root(z^3 - m*z^2 - n/24, z, 3)

       可见此方程并不存在解析解，即未知数x并不能写成关于m和n的表达式。root为matlab中自带的求根函数，说明在已知m，n的情况下，x的解即为z^3 - m*z^2 - n/24关于z的解，由于三次方程其解必然为三个，故需在m和n都已知的情况下，才能计算出x的数值解。那么对于此公式，写成函数形式的代码时可以直接将运行结果的语句复制过去，再使用double函数。

代码如下：

function cgm2 =  Status(E,r1,r2,t1,t2,alpha,cgm1)
clear all；
syms z;
l=50:50:700;

for i = 1:14
    
    m=cgm1-E*r1^2*l(i)^2/(24*cgm1^2)-alpha*E*(t2-t1);
    n=E*r2^2*l(i)^2;
    
    slove=root(z^3 - m*z^2 - n/24, z, 1); 
    cgm2(i)=double(slove);
    
    
end

       此论文中需要得到L从50到700每次递加50的每个情况下的σ02，故写成这个形式。

命令行中调用：cgm2=Status(65000,31.13e-3,31.13e-3,-20,40,20.5e-6,91.56)

运行结果：

cgm2 =

  1 至 6 列

   23.1115   32.8301   40.5246   46.8721   52.2131   56.7579

  7 至 12 列

   60.6533   64.0097   66.9137   69.4353   71.6323   73.5525

  13 至 14 列

   75.2360   76.7166

       按照这个方法，可以解决一般论文中的公式求解。再比如下面的例子：

θT-1eθT=θ2λ(h+θC)-1

       除了T以外其余参数可知，同样的，将等号右边整合为m，θT等于n，尝试找n的解析式，在最后写函数时再求出n再除去θ。

代码：

syms m n;

f1 = [ (n-1)*exp(n)==m];

solve(f1,n)

运行结果：lambertw(0, m*exp(-1)) + 1

       lambertw也是matlab中自带的函数，W = lambertw(X)其意为：We^w=X。事实上，我们此时无须具体的求解出n，只需得到此次运行结果，之后直接复制到函数中即可。

代码如下：

function T =  Status(ct,K,lamd,C,h)
clear all；
syms m n;

m = ct^2*K/(lamd*(h+ct*C))-1;
f1 = [ (n-1)*exp(n)==m];

ans=solve(f1,n);

T=double(ans)/ct;
  
end

在命令行调用：T=Status(0.05,300,2000,10,2)

运行结果：

T =

    0.3444

       可见用这种方法可以很简单的实现论文中的公式，首先尝试获取在matlab中关于公式中某一变量的表达式，再使用double函数取得其数值解，即可写成某一变量关于其余变量的函数。