今年连续更新了4篇水博客。现在我都不敢在开口说什么新的新年计划了,之前的都只完成了一半。自控能力要加强,不要耽误时间。虽然嘴上说着996不好,但是不抓紧时间,吃亏的只能是自己。去年没做完的整理上学期间学习的代码工作在本年内会继续。之前该补的东西,尽量都补上,要不都不好意思开新的坑。

除此之外,还有一条新的经验要牢记:永远不到等到文件夹乱得找不到文件了,再去整理,养成良好的电脑文件归类习惯很重要。先乱放,再整理,是在是太耽误时间了。

人的一生是熵不断增加的过程。整理好自己的用品,安排好各种时间,能降低自己平时生活中的熵(混乱程度)。

再多写一些展望吧。。

  • 我最喜欢的 RNA 将在医疗领域发挥更重大的作用(不止是疫苗)。
  • 如何让其他星球产生含氧且适宜人类生存的大气层是个研究方向,感觉还没准能成功呢。
  • 人机接口这东西,最近还最好先关注一些肢体神经,饶了大脑吧。

互联网公司职工按 996 甚至每月休 2 天的频率来开展工作,因身体和精神原因出事频繁,尤其某个“x多多”近期频繁上新闻。 最近读了马克思的《1844年经济学哲学手稿》,认识到某地还是个“初级”的“资本主义”社会。简直是天天抓理论学习,但永远都是学一套、说一套、做一套。

在《1844年经济学哲学手稿》中马克思提出了异化劳动的观点。异化的最初含义可以理解为“自我外化为非我,从而使原来与自我同一的东西,变成异己的东西”。例如,我创作一件木雕作品,把自己的思想和手艺固化到这个作品里,也就是将自己的“部分思想”和“手艺”同自己分离,形成一件雕塑,变成不是自己的的东西。马克思的观点认为:异化作为社会现象,是与阶级一起产生的,是人的物质生产机器产品变成异己力量,反过来统治人的社会现象。私有制是异化产生的社会根源,社会分工固定是它的自然缘由。在异化活动中,人的能动性丧失了。遭到异己的物质力量或精神力量的奴役,从而使人的个性不能全面发展,只能片面发展,甚至畸形发展。在资本主义社会,异化达到最严重的程度。

资本主义的工厂中,让工人延长上班时间,生产更多的产品,产生更多的利润。但产品的利润归资本家所有,反过来工人只能通过自己劳动所的工资来购买产品。工人为了购买必需品,而进行劳动,这时候,劳动就不是工人“自由自觉”的活动,具有强制性、外在性和异己性。这同互联网行业的996,以及某些实验室的研究生培养,没有任何区别。互联网行业里,大家每天回家就睡个觉,其余时间都在通勤、搬砖、吃喝拉撒,自己的产品获得的利润归资本家所有,而自己只能通过获得的工资来购买产品,期权这些等着解禁还不知要熬过多少年,前提是身体健康情况允许你继续熬。这也是互联网行业35岁以上不要的原因,当不了电池继续发电,自然要淘汰。一些打卡严格的高等教育实验室也是严格规定学生的工作时间,本来收入就只能够吃饭,还可能时常受到精神打击。在这些异化活动里,人的能动性完全丧失了,绝大部分人只能片面的发展。纵观世界,东亚三国,两个资本主义国家、一个社会主义国家,都差不多这样。还真说不好这么“卷”是因为我们自身的发达程度不够,还是人口原因或者文化原因。

马克思根据唯物主义辩证法,唯物史观和剩余价值学说,还论证过资本主义必然灭亡,我不知他看到现今的情况,会不会发表什么新的看法,毕竟他老家的工人待遇比这边的好太多。在发达国家里比较的话,估计 Musk 的公司要排除在外,我觉得他家公司的工作时长应该和我们的 996 有一拼。


一阶常微分方程初值问题,一般表示如下所示:

\[\begin{cases} y^{'}=f(x,y)\\ y(x_{0})=y_{0} \end{cases}\]

欧拉法,改进的欧拉法,龙格-库塔法都是基于同样的原理,即用切线\(y^{'}=f(x,y)\)去逼近原方程的曲线\(y=f(x)\)1

Euler_Method

那么怎么作出切线呢,\(y^{'}=f(x,y)\)就是这个切线的方程。

欧拉法计算时,令每次前进的步长为\(h\),每次根据步长h根据切线\(f(x,y)\)的方向求下一个\(y\)的位置。

\[\begin{cases} y_{i+1}=y_{i}+hK_{1}\\ K_{1}=f(x_{i},y_{i}) \end{cases}\]

改进的欧拉法,就是在欧拉法的基础上修改方向(斜率),是得每步计算的切线更贴近原曲线。

\[\begin{cases} y_{i+1}=y_{i}+\frac{h}{2}(K_{1}+K_{2})\\ K_{1}=f(x_{i},y_{i})\\ K_{2}=f(x_{i+1},y_{i}+hK_{1}) \end{cases}\]

如果设法在\([x_{i},x_{i+1}]\)内多预报几个点的斜率值,然后将它们加权平均作为平均斜率\(K^{*}\),则有可能构造出更高精度的计算格式,龙格-库塔法就是在上述改进欧拉法的基础上,继续构造新的\(K\),来达到更加精确的逼近原曲线的目的。

一般常用的方法是四阶龙格-库塔法,计算公式如下:

\[\begin{cases} K_{1}=f(x,y)\\ K_{2}=f(x_{i+\frac{1}{2}},y_{i}+\frac{h}{2}K_{1})\\ K_{3}=f(x_{i+\frac{1}{2}},y_{i}+\frac{h}{2}K_{2})\\ K_{4}=f(x_{i+1},y_{i}+hK_{3})\\ y_{i+1}=y_{i}+\frac{h}{6}(K_{1}+2K_{2}+2K_{3}+K_{4}) \end{cases}\]

四阶龙格-库塔法程序如下所示。

#include"stdio.h"
#include"stdlib.h"
void fun1(double a,double b,double h)
{  double x,y=0,i,k1,k2,k3,k4;
   for(i=1;i<=((b-a)/h)+1;i++)
 { x=i*h;
   k1=h*(1-y);
   k2=h*(1-(y+1/2.0*k1));
   k3=h*(1-(y+1/2.0*k2));
   k4=h*(1-(y+k3));
   y=y+1/6.0*(k1+2*k2+2*k3+k4);
 }  
 printf("%lf\n",y);
      }
      
void fun2(double a,double b,double h)
{  double x,y=1,i,k1,k2,k3,k4;   
   for(i=0;i<=((b-a)/h);i++)
 { x=i*h;
   k1=h*(x*y*y);
   k2=h*((x+h/2)*(y+1/2.0*k1)*(y+1/2.0*k1));
   k3=h*((x+h/2)*(y+1/2.0*k2)*(y+1/2.0*k2));
   k4=h*((x+h)*(y+k3)*(y+k3));
   y=y+1/6.0*(k1+2*k2+2*k3+k4);
 }  
 printf("%lf\n",y);
      }

int main()
{double a=0,b=1,h=0.1;
      fun1(a,b,h);
      fun2(a,b,h);
      }

1. 雅克比(Jacobi)迭代法求方程组的根

雅克比迭代法非常简单,对于一个给定的n×n方程组\(\displaystyle A\mathbf {x} =\mathbf {b}\)

\[\displaystyle A={\begin{bmatrix}a_{11}a_{12}\cdots a_{1n}\\a_{21}a_{22}\cdots a_{2n}\\\vdots \vdots \ddots \vdots \\a_{n1}a_{n2}\cdots a_{nn}\end{bmatrix}},\qquad \mathbf {x} ={\begin{bmatrix}x_{1}\\x_{2}\\\vdots \\x_{n}\end{bmatrix}},\qquad \mathbf {b} ={\begin{bmatrix}b_{1}\\b_{2}\\\vdots \\b_{n}\end{bmatrix}}.\]

可以把 A分解成两个矩阵

\[\displaystyle A=D+R\qquad \qquad D={\begin{bmatrix}a_{11}0\cdots 0\\0a_{22}\cdots 0\\\vdots \vdots \ddots \vdots \\00\cdots a_{nn}\end{bmatrix}},\qquad R={\begin{bmatrix}0a_{12}\cdots a_{1n}\\a_{21}0\cdots a_{2n}\\\vdots \vdots \ddots \vdots \\a_{n1}a_{n2}\cdots 0\end{bmatrix}}\]

方程组可以改写为 \(\displaystyle D\mathbf {x} =\mathbf {b} -R\mathbf {x}\),该公式即为迭代公式。

又可写成 \(\displaystyle \mathbf {x} ^{(k+1)}=D^{-1}(\mathbf {b} -R\mathbf {x} ^{(k)})\)。

注意这种解方程组的问题,矩阵都要满一定的收敛条件,在这里是迭代矩阵的谱半径小于1,\(\rho (D^{-1}R)<1\)。

#include"math.h"
#include"stdio.h"
#include"stdlib.h"
void main()
{
	double x1,x2,x3,x11,x22,x33;
    int i=0;
    printf("No.0008雅克比迭代法求方程组的根\n\n输入初值:");
    scanf("%lf %lf %lf",&x1,&x2,&x3);
    while(fabs(x11-(1-x2-x3)/(-8.0))>0.001 && fabs(x22-(16-x1-x3)/(-5.0))>0.001 && fabs(x33-(7-x1-x2)/(-4.0))>0.001)
  {
	x11=(1-x2-x3)/(-8.0);//写迭代式 
	x22=(16-x1-x3)/(-5.0);//写迭代式
	x33=(7-x1-x2)/(-4.0);//写迭代式
	x1=x11;
	x2=x22;
	x3=x33;
    i++;
  }	
	printf("解为:\nx1=%lf\nx2=%lf \nx3=%lf\n迭代次数:%d\n",x1,x2,x3,i);
system("pause");
} 

2. 高斯-赛德尔(Gauss-Seidel)迭代法求方程组的根

高斯-赛德尔迭代法和雅克比迭代法的区别在于,将A分解为一个上三角矩阵和下三角矩阵的形式。

\[\displaystyle A=L_{*}+U\qquad {\text{其中}}\qquad L_{*}={\begin{bmatrix}a_{11}0\cdots 0\\a_{21}a_{22}\cdots 0\\\vdots \vdots \ddots \vdots \\a_{n1}a_{n2}\cdots a_{nn}\end{bmatrix}},\quad U={\begin{bmatrix}0a_{12}\cdots a_{1n}\\00\cdots a_{2n}\\\vdots \vdots \ddots \vdots \\00\cdots 0\end{bmatrix}}\]

这样,代入迭代式后为:

\[\displaystyle L_{*}\mathbf {x} =\mathbf {b} -U\mathbf {x}\] \[\displaystyle \mathbf {x} ^{(k+1)}=L_{*}^{-1}(\mathbf {b} -U\mathbf {x} ^{(k)})\]

系数矩阵 A 严格对角占优或对称正定时,高斯-赛德尔迭代法必收敛。

#include"math.h"
#include"stdio.h"
#include"stdlib.h"
main()
{
	double x1,x2,x3,x11,x22,x33;
    int i=0;
    printf("No.0009高斯赛德尔(Gauss-Seidel)迭代法求方程组的根\n\n输入初值:");
    scanf("%lf %lf %lf",&x1,&x2,&x3);
    while(fabs(x11-(1-x2-x3)/(-8.0))>0.01 && fabs(x22-(16-x1-x3)/(-5.0))>0.01 && fabs(x33-(7-x1-x2)/(-4.0))>0.01)
  {
	x11=x1;x22=x2;x33=x3;
	x1=(1-x2-x3)/(-8.0);//写迭代式 
	x2=(16-x1-x3)/(-5.0);//写迭代式
	x3=(7-x1-x2)/(-4.0);//写迭代式
    i++;
  }	
	printf("解为:\nx1=%lf\nx2=%lf \nx3=%lf\n迭代次数:%d\n",x1,x2,x3,i);
system("pause");
}

3. 牛顿迭代法求方程组的根

这个根牛顿迭代法求解方程的根一样。。。

\[\displaystyle x_{n+1}=x_{n}-J_{F}(x_{n})^{-1}F(x_{n})\]

其中 \(J_{F}(x)\)是雅克比矩阵。

#include"math.h"
#include"stdio.h"
#include"stdlib.h"
main()
{
	double x1,x2,x11,x22;
    int i=0;
    printf("No.0010牛顿迭代法求方程组的根\n\n输入初值:");
    scanf("%lf %lf",&x1,&x2);
	x11=x1-(2*x2/(2*x2-8*x1)*(x1+2*x2-3)-2/(2*x2-8*x1)*(2*x1*x1+x2*x2-5));//写迭代式 
	x22=x2-(-4*x1/(2*x2-8*x1)*(x1+2*x2-3)+1/(2*x2-8*x1)*(2*x1*x1+x2*x2-5));
    i++;//写迭代式
 while(fabs(x11-x1)>0.001 && fabs(x22-x2)>0.001)
  { x1=x11;
	x2=x22;
	x11=x1-(2*x2/(2*x2-8*x1)*(x1+2*x2-3)-2/(2*x2-8*x1)*(2*x1*x1+x2*x2-5));//写迭代式 
	x22=x2-(-4*x1/(2*x2-8*x1)*(x1+2*x2-3)+1/(2*x2-8*x1)*(2*x1*x1+x2*x2-5));//写迭代式
    i++;
  }	
	printf("解为:\nx1=%lf\nx2=%lf \n迭代次数:%d\n",x11,x22,i);
system("pause");
} 

1. 二分法求方程的根

二分法可行的条件是有介值定理支持。

介值定理说明在实数范围内,\([a,b]\)区间上可以画出一个连续曲线。\(f(a)<u<f(b)\) 存在,那么存在\(c\), \(a<c<b\),使得\(f(c)=u\)。

二分法就是\(a\)与\(b\)异号时,介值定理的特殊情况,即存在实数\(c\)使得\(f(c)=0\)存在。

#include"math.h"
#include"stdio.h"
#include"stdlib.h"
double fun(double k)
{return 1-k-sin(k);//所求方程 
}
void main()
{double a,b,x;int i=0;
printf("No.0004二分法求方程的根\n\n输入区间左端点:");
scanf("%lf",&a);
printf("输入区间右端点:"); 
scanf("%lf",&b);
x=(a+b)/2.0;
while(/*fabs(fun(x))>0.00005||*/fabs(a-b)>0.00005)//近似范围 
{
if(fun(a)*fun(x)<0) b=x;
if(fun(b)*fun(x)<0) a=x;
x=(a+b)/2.0;
i++;}
printf("方程的根:%lf\n",x);
printf("迭代次数:%d\n\n",i);
system("pause");
}

2. 迭代法求方程的根

迭代法主要是将\(f(x)=0\)这种等式,该写成\(x=\varphi(x)\),然后给定一个初始的\(x_0\),根据公式\(x_1=\varphi(x_0)\),求出\(x_1\),接下来,比较\(x_1\)和\(x_0\)的差值大小,如果差值等于0或者一个非常非常小的数,那么则求出最终的\(x\)即为\(x_1\);如果差值还很大,则进行迭代步骤,将\(x_1\)带入到公式右侧\(\varphi(x_1)\),求出左侧的\(x_2\),重复上述比较步骤,直到\(x\)收敛于一个定值。

#include"math.h"
#include"stdio.h"
#include"stdlib.h"
double fun(double k)
{
	return pow(1+k*k,1/3.0);//所求方程 
}
void main()
{double a,b,x,x1,l=1;
int i=0,p;
printf("No.0005迭代法方程的根\n\n选择输入数据种类:\n1.数值\n2.区间\n");
printf("\n");
while(l==1)
{scanf("%d",&p);
if(p!=1&&p!=2) printf("请输入指定序号\n\n");//判断输入的序号是否正确 
else l=0;
}
switch(p)
 {
 case 1:
		printf("数据种类\n1.数值\n输入数值:");
		scanf("%lf",&x);
		x1=fun(x);i=i+1;
		break;
 case 2:printf("数据种类\n2.区间\n");
		printf("输入区间左端点:");
		scanf("%lf",&x);
		printf("输入区间右端点:");
		scanf("%lf",&x1);x1=fun((x+x1)/2.);i=i+1;
		break;
 }
while(fabs(x-x1)>0.00005)
{
x=fun(x1);
i=i+1;
x1=fun(x);
i=i+1;
}
printf("方程的根:%lf\n",x1); 
printf("迭代次数:%d\n\n",i);
system("pause");
}

3.牛顿迭代法

牛顿迭代法又称为牛顿-拉弗森方法(Newton-Raphson method),它比一般的迭代法有更高的收敛速度。

Newton-Raphson method

\(f(x_0)\)点的切线是\(f(x)\)的线性逼近。离\(f(x_0)\)点距离越近,这种逼近的效果也就越好,也就是说,切线与曲线之间的误差越小。所以我们可以说在\(f(x_0)\)点附近,\(切线\approx f(x)\)。

牛顿-拉弗森方法提出来的思路就是利用切线是曲线的线性逼近这个思想,随便找一个曲线上的\(f(x_0)\)点(为什么随便找,根据切线是切点附近的曲线的近似,应该在根点附近找,但是很显然我们现在还不知道根点在哪里),做一个切线,切线的根\(x_1\)(就是和x轴的交点)与曲线的根,还有一定的距离。我们从这个切线的根\(x_1\)出发,做一根垂线,和曲线相交于\(f(x_1)\)点,继续重复刚才的工作,多次迭代后会越来越接近曲线的根,即迭代收敛。

在计算时,首先选择一个\(x_0\),然后计算相应的\(f(x_0)\)和切线斜率\(f^{'}(x_0)\),接下来计算穿过点(\(x_0\),\(f(x_0)\))并且斜率为\(f^{'}(x_0)\)的直线和\(x\)轴的交点的\(x_1\)坐标,公式为: \(0=(x_1-x_0)f^{'}(x_0)+f(x_0)\), 将其转化为\(x_1=\varphi(x_0)\)的形式,即有:\(x_{n+1}=x_n-\frac{f(x_n)}{f^{'}(x_n)}\)

将初始的\(x_0\)代入到上述公式中,经过多次迭代,可以求得最后的根。

#include"math.h"
#include"stdio.h"
#include"stdlib.h"
double fun(double k)
{
	return k*k*k-k*k-1;//所求方程 
}
double fun1(double k)
{  
	return 3*k*k-2*k;//所求方程的一阶导函数 
}
void main()
{double a,b,x,x1,l=1;
int i=0,p;
printf("No.0006牛顿迭代法求方程的根\n\n选择输入数据种类:\n1.数值\n2.区间\n");
printf("\n");
while(l==1)
{scanf("%d",&p);
if(p!=1&&p!=2) printf("请输入指定序号\n\n");//判断输入的序号是否正确 
else l=0;
}
switch(p)
{
 case 1:
		printf("数据种类\n1.数值\n输入数值:");
		scanf("%lf",&x1);
		break;
 case 2:printf("数据种类\n2.区间\n");
		printf("输入区间左端点:");
		scanf("%lf",&x);
		printf("输入区间右端点:");
		scanf("%lf",&x1);x1=(x+x1)/2.;
		break;
 }
while(fabs(fun(x1)/fun1(x1))>0.00005)
{
x=x1-fun(x1)/fun1(x1);
x1=x;
i++;
}
printf("方程的根:%lf\n",x1); 
printf("迭代次数:%d\n\n",i);
system("pause");
}

4. 弦割法

弦割法(Secant Method)是基于牛顿法的一种改进,基本思想是用弦的斜率近似代替目标函数的切线斜率,并用割线与横轴交点的横坐标作为方程式的根的近似。

SecantMethod

公式就不推了,参见维基百科或者百度百科。

#include"math.h"
#include"stdio.h"
#include"stdlib.h"
double fun(double k)
{
	return k*k*k-k*k-1;//所求方程 
}
main()
{double a,b,x,x1,l=1,x0,x2;
int i=0,p;
printf("No.0007弦割法求方程的根\n\n选择输入种类:\n1.单点\n2.双点\n");
printf("\n");
while(l==1)
{scanf("%d",&p);
if(p!=1&&p!=2) printf("请输入指定序号\n\n");//判断输入的序号是否正确 
else l=0;
}
switch(p)
{
 case 1:
		printf("种类\n1.单点\n");
		printf("输入点1:");
		scanf("%lf",&x);
		printf("输入点2:");
		scanf("%lf",&x1);
        while(fabs(fun(x1)*(x1-x)/(fun(x1)-fun(x)))>0.00005)
        {
        x2=x1-fun(x1)*(x1-x)/(fun(x1)-fun(x));
        x1=x2;
        i++;
        }
        printf("方程的根:%lf\n",x1); 
        printf("迭代次数:%d\n\n",i);
		break;
 case 2:printf("种类\n2.双点\n");
		printf("输入点1:");
		scanf("%lf",&x);
		printf("输入点2:");
		scanf("%lf",&x1);
        while(fabs(fun(x1)*(x1-x)/(fun(x1)-fun(x)))>0.00005)
        {
        x2=x1-fun(x1)*(x1-x)/(fun(x1)-fun(x));
        x=x1;
        x1=x2;
        i++;
        }
        printf("方程的根:%lf\n",x1); 
        printf("迭代次数:%d\n\n",i);
	    break;
}
system("pause");