正者，正也。其心以为不然者，天门弗开矣。作为《深度学习入门》的阅读笔记，本文简略概述深度学习的学习过程，损失函数，梯度法相关的内容，并辅以python实现需要事先了解的本文并不是从零开始的，需要实现了解蜻蜓点水python蜻蜓点水python_dlc简略感知机简易神经网络_推理和正向传播不同的学习首先，学习有不同的类型，人想，机器学习，深度学习；然后，学习的基础是特征向量，机器学习中特征量还是人想出来的，深度学习则没有这一步骤，基本实现了全程机器管理损失函数有道是能观测就能干涉，能干涉就能控制；那么应该观测什么变量才能知道学习的进度，从而干涉和控制学习的过程呢？这里直接给出答案：损失函数。损失函数理论上可以使用任何函数，这里使用均方误差和交叉熵误差 。这里有一个问题，为什么不使用精度（即准确度）来最为指标呢？首先说结论：在进行神经网络的学习时，不能将识别精度作为指标，并且对于大部分导数（后面会提到）为0的函数，都不能作为指标，会导致学习无法进行（无法更新参数），原因在了解学习的过程后自然能理解。均方误差直接公式，实现如下：y_k是前向传播的输出t_k是监督数据的标签（一般是one-hot表示）k是维度import numpy as np def mean_squared_error(y, t): return 0.5 * np.sum((y-t)**2) print(mean_squared_error(np.array([0.1,0.7,0.2]),np.array([0,1,0]))) #判断正确 print(mean_squared_error(np.array([0.1,0.7,0.2]),np.array([1,0,0]))) #判断错误 #0.07000000000000002 #0.67可以看出，如果判断正确，损失函数会是一个较小的数字，如果判断失败，损失函数会是一个较大的数字。交叉熵误差直接公式，实现如下：y_k是前向传播的输出t_k是监督数据的标签（一般是one-hot表示）k是维度def cross_entropy_error(y, t): delta = 1e-7 return -np.sum(t * np.log(y + delta)) print(cross_entropy_error(np.array([0.1,0.7,0.2]),np.array([0,1,0]))) print(cross_entropy_error(np.array([0.1,0.7,0.2]),np.array([1,0,0]))) #0.3566748010815999 #2.302584092994546其中delta是为了避免log的参数有0的情况Mini-Batch学习和推理的过程一样，学习也可以用批处理的方法完成，不过实现上有些区别；这里以交叉熵误差为例第一个是损失函数公式的变化：可见其实就是取平均，实现如下：def cross_entropy_error2(y, t): batch_size = 1 if y.ndim > 1: batch_size = y.shape[0] return -np.sum(t * np.log(y + 1e-7)) / batch_size #测试不同维度 print(cross_entropy_error2(np.array([0.1,0.7,0.2]),np.array([0,1,0]))) print(cross_entropy_error2(np.array([0.1,0.7,0.2]),np.array([1,0,0]))) print(cross_entropy_error2(np.array([[0.1,0.7,0.2],[0.1,0.7,0.2]]),np.array([[0,1,0],[0,1,0]]))) print(cross_entropy_error2(np.array([[0.1,0.7,0.2],[0.1,0.7,0.2]]),np.array([[1,0,0],[1,0,0]])))第二个是训练数据其实还是很多的，想一次性全部学习是高耗能甚至有些不现实的，所以需要使用np.random.choice抽取一个batch来学习# 这里的路径按本地实际来 from dfs.dataset.mnist import load_mnist import sys,os sys.path.append(os.pardir) (x_train, t_train), (x_test, t_test) = load_mnist(flatten=True, normalize=False,one_hot_label=True) # 输出各个数据的形状 print(x_train.shape) # (60000, 784) # 输出各个数据的形状 print(t_train.shape) # (60000,) print(x_test.shape) # (10000, 784) print(t_test.shape) # (10000,) train_size = x_train.shape[0] batch_size = 10 batch_mask = np.random.choice(train_size, batch_size) print(batch_mask) # 输出抽选的索引 x_batch = x_train[batch_mask] t_batch = t_train[batch_mask] print(x_batch.shape) # (10, 784) print(t_batch.shape) # (10,) #(60000, 784) #(60000, 10) #(10000, 784) #(10000, 10) #[ 6681 4809 12952 4669 1444 20997 56685 39195 24206 8649] #(10, 784) #(10, 10)第三如果标签不是one-hot表示，有那么一个特殊的实现方法，思路是虽然是求和，但只有正确项会有值def cross_entropy_error2(y, t): if y.ndim == 1: t = t.reshape(1, t.size) y = y.reshape(1, y.size) batch_size = y.shape[0] print(y[np.arange(batch_size), t]) return -np.sum(np.log(y[np.arange(batch_size), t] + 1e-7)) / batch_size print(cross_entropy_error2(np.array([[0.1,0.7,0.2],[0.1,0.7,0.2]]),np.array([1,2]))) #[0.7 0.2] #0.9830561067579126另外有个概念会在后面用到这里先提出：epoch是一个单位。一个epoch表示学习中所有训练数据均被使用过一次时的更新次数。比如，对于 10000笔训练数据，用大小为 100笔数据的mini-batch进行学习时，重复随机梯度下降法100次，所有的训练数据就都被“看过”了 。此时，100次就是一个epoch。导数和梯度有了损失函数，就可以观测学习的进度，但是如何在此基础上对其进行控制呢？要做到这点需要解决三个问题，函数的值是会通过自变量的变化而变化的，改哪个变量会使损失函数变小？如何改？要改多少？直接说结论：修改的变量是权重和偏置，依靠导数计算出的梯度来修改自变量，改多少依赖于人工设定的学习率。导数不准确的说，导数就是x周围很小范围的斜率（中心差分），既不是x左边的斜率也不是x右边的斜率（前向差分）；数值微分和真实的解析解必然是有所差距的。求导数的实现如下：def numerical_diff(f,x,h=1e-4): return (f(x+h)-f(x-h))/(2*h) def func_test(x): return 2*x+1 numerical_diff(func_test,np.array([2,3]))偏导数和梯度多个变量的函数的导数就是偏导数，当然一个偏导数只能对其中一个变量来做，另外的变量就当作常数对待，称为函数对该变量的偏导数，这是解析解的做法；对于数值解则更加的简单一点，固定其他变量，只变化一个变量来求导数就是函数针对这个变量的偏导数。把所有变量的偏导数组成向量，就是梯度，梯度指向的方向是函数变化（增加）率最大的方向。所以对其取负数，就是损失函数下降最快的方向。实现如下：def numerical_gradient(f, x): h = 1e-4 # 0.0001 grad = np.zeros_like(x) it = np.nditer(x, flags=['multi_index'], op_flags=['readwrite']) while not it.finished: idx = it.multi_index tmp_val = x[idx] x[idx] = float(tmp_val) + h fxh1 = f(x) # f(x+h) x[idx] = tmp_val - h fxh2 = f(x) # f(x-h) grad[idx] = (fxh1 - fxh2) / (2*h) x[idx] = tmp_val # 还原值 it.iternext() return grad def test_func(x): return x[0]*x[0]+x[1]*x[1] numerical_gradient(test_func,np.array([3.0,4.0])) #array([6., 8.])NumPy 迭代器对象 numpy.nditer 提供了一种灵活访问一个或者多个数组元素的方式。见NumPy 迭代数组 | 菜鸟教程梯度下降法按上述说的，按照梯度反方向更新参数就是梯度下降法，公式如下：其中偏导数前面的常数就是学习率，人为指定，规定学习的程度，过大过小都不好，这里默认0.01，实现如下：def gradient_descent(f,init_x,lr=0.01,step_num=100): x = init_x for i in range(step_num): grad = numerical_gradient(f,x) x -= lr*grad return x def function_2(x): return x[0]*x[0]+x[1]*x[1] init_x=np.array([-3.0,4.0]) gradient_descent(function_2,init_x,lr=0.1,step_num=100) #array([-6.11110793e-10, 8.14814391e-10])另外，像这种人为指定的参数，称为超参数，像是mini-batch的大小，输入的形状，学习的次数都是超参数总结和实现现在就可以继续上章，完成MINIST的学习步骤了，但是还有一些细节需要说明神经网络的梯度法在实际的神经网络中要求的是损失函数对于权重的梯度，也就是对权重矩阵中的每个成员求梯度，最后结果还是矩阵：对于实现而言，只要能保证矩阵的形状没有问题即可学习算法的步骤步骤4则是重复1-3的步骤而这种随机选取+计算梯度的方法称为：随机梯度下降法（SGD）实现先列出代码# 这里的路径按本地实际来 import numpy as np from dfs.dataset.mnist import load_mnist import sys,os sys.path.append(os.pardir) (x_train, t_train), (x_test, t_test) = load_mnist(flatten=True, normalize=False,one_hot_label=True) # 输出各个数据的形状 #print(x_train.shape) # (60000, 784) # 输出各个数据的形状 #print(t_train.shape) # (60000,10) #print(x_test.shape) # (10000, 784) #print(t_test.shape) # (10000,10) train_size = x_train.shape[0] batch_size = 10 batch_mask = np.random.choice(train_size, batch_size) #print(batch_mask) x_batch = x_train[batch_mask] t_batch = t_train[batch_mask] print(x_batch.shape) # (10, 784) print(t_batch.shape) # (10,10) def cross_entropy_error(y, t): batch_size = 1 if y.ndim > 1: batch_size = y.shape[0] return -np.sum(t * np.log(y + 1e-7)) / batch_size def softmax(a): c = np.max(a) exp_a = np.exp(a-c) #溢出对策 sum_exp_a = np.sum(exp_a) y = exp_a /sum_exp_a return y def sigmoid(x): return 1/(1+np.exp(-x)) class SigMoidLayer: def __init__(self,w,b,s="SigMoid"): self.layer_name = s self.w=w self.b = b def forward(self,x): a=np.dot(x,self.w)+self.b return sigmoid(a) class OutputLayer2: def __init__(self,w,b,s="Ouput"): self.layer_name = s self.w=w self.b = b def forward(self,x): a=np.dot(x,self.w)+self.b return softmax(a) def numerical_gradient(f, x): h = 1e-4 # 0.0001 grad = np.zeros_like(x) it = np.nditer(x, flags=['multi_index'], op_flags=['readwrite']) while not it.finished: idx = it.multi_index tmp_val = x[idx] x[idx] = float(tmp_val) + h fxh1 = f(x) # f(x+h) x[idx] = tmp_val - h fxh2 = f(x) # f(x-h) grad[idx] = (fxh1 - fxh2) / (2*h) x[idx] = tmp_val # 还原值 it.iternext() return grad class TestMNIST: # 正常情况下，权重或者偏置都应该是学习的结果，这里直接赋值只是一种假设 def __init__(self): self.rate=0.01 self.grads_b = {} self.grads_w = {} self.network_w = {} self.network_b = {} self.network_w['L1'] = np.random.randn(784,50) self.network_b['L1'] = np.zeros(50) self.network_w['L2'] = np.random.randn(50,100) self.network_b['L2'] = np.zeros(100) # 输出层 self.network_w['L3'] = np.random.randn(100,10) self.network_b['L3'] = np.zeros(10) self.layers=['L1','L2'] def forward(self,x): tmp=x for l in self.layers: tmp = SigMoidLayer(self.network_w[l],self.network_b[l],l).forward(tmp) #print(l,tmp) tmp = OutputLayer2(self.network_w['L3'],self.network_b['L3'],'L3').forward(tmp) return tmp def loss(self,x,t): return cross_entropy_error(self.forward(x),t) def grad(self,x,t): loss_W = lambda W: self.loss(x, t) self.grads_w={} self.grads_b={} self.grads_w['L1']=numerical_gradient(loss_W,self.network_w['L1']) self.grads_b['L1']=numerical_gradient(loss_W,self.network_b['L1']) print('L1 Done.') self.grads_w['L2']=numerical_gradient(loss_W,self.network_w['L2']) self.grads_b['L2']=numerical_gradient(loss_W,self.network_b['L2']) print('L2 Done.') self.grads_w['L3']=numerical_gradient(loss_W,self.network_w['L3']) self.grads_b['L3']=numerical_gradient(loss_W,self.network_b['L3']) print('L3 Done.') def update(self): self.network_w['L1']-=self.rate*self.grads_w['L1'] self.network_b['L1']-=self.rate*self.grads_b['L1'] self.network_w['L2']-=self.rate*self.grads_w['L2'] self.network_b['L2']-=self.rate*self.grads_b['L2'] self.network_w['L3']-=self.rate*self.grads_w['L3'] self.network_b['L3']-=self.rate*self.grads_b['L3'] nett = TestMNIST() #测试代码 #nett.forward(np.random.randn(100,784)) # 学习过程 mini_batch_size = x_train.shape[0] select_size = 100 loss_point =[] train_time=10000 for i in range(train_time): # 选择数据 mask = np.random.choice(mini_batch_size, select_size) x_batch = x_train[mask] t_batch = t_train[mask] print(x_batch.shape,t_batch.shape) print("Loop",i+1) # 损失函数 loss_p=nett.loss(x_batch,t_batch) loss_point.append(loss_p) print("Loss:",loss_p) # 计算梯度 nett.grad(x_batch,t_batch) # 更新参数 nett.update() #推理过程 batch_size=100 accuracy_cnt=0 for i in range(0,len(x_test),batch_size): x_batch = x_test[i:i+batch_size] t_batch = t_test[i:i+batch_size] y_batch = nett.forward(x_batch) #print(y_batch.shape) p=np.argmax(y_batch,axis=1) t=np.argmax(t_batch,axis=1) accuracy_cnt+=(np.sum(p==t)/batch_size) print("Accuracy:" + str(float(accuracy_cnt) / (len(x_test)/batch_size))) import matplotlib.pyplot as plt # 生成数据 xx = np.arange(0, train_time) # 绘制图形 plt.plot(xx, loss_point) plt.show()首先要注意的是矩阵的形状一定要“严丝合缝”，这是做矩阵计算的基础np.zeros(50)是创建0矩阵的方法，第一个参数是shape然后当运行该代码时可以发现一点，就是随机梯度下降的处理速度十分缓慢，所以在之后会需要反向传播这一方法来加速梯度的计算评价方法评价一个神经网络的方法就是评价它的泛化能力，也就是对于陌生数据的识别精度：如果只对训练数据识别精度较好，则是发生了过拟合一般来说，当学习经过一个epoch，就可以对测试数据和训练数据进行推理，然后进行可视化对比，这里只给出伪代码：epoch是一个单位。一个epoch表示学习中所有训练数据均被使用过一次时的更新次数。比如，对于 10000笔训练数据，用大小为 100笔数据的mini-batch进行学习时，重复随机梯度下降法100次，所有的训练数据就都被“看过”了。此时，100次就是一个epoch。不幸的结果和原著对比，这里的效果差了很多（最重要的是运行时间巨长！！！！整整2个星期），原因会放在未来探索，可能的原因一方面是初始化参数的问题，一方面可能是层级结构的问题。准确度在：Accuracy:0.7011999999999997 而不是0.9以上，梯度下降的也不太对，如下图。但不管怎么样，SGD是有效果的，接下来需要解决其速度太慢的问题。
化其万物而不知其禅之者，焉知其所终？焉知其所始？正而待之而已耳。作为《深度学习入门》的阅读笔记，本文简略概述正向传播和激活函数的内容，并辅以python实现需要事先了解的本文并不是从零开始的，需要实现了解蜻蜓点水python蜻蜓点水python_dlc简略感知机学习和推理和求解机器学习问题的步骤（分成学习和推理两个阶段进行）一样，使用神经网络解决问题时，也需要首先使用训练数据（学习数据）进行权重参数的学习；进行推理时，使用刚才学习到的参数，对输入数据进行分类。 首先，通过机器学习解决问题通常分为两个阶段：学习和推理；对于感知机或者神经网络，学习的过程就是调整参数的过程，是需要事先执行的步骤；在完成学习后，使用确定的参数解决实际问题的过程，称为推理；在神经网络或感知机中，通过前向传播来实现；从感知机到神经网络回顾下异或问题，在感知机的实现中，输入先通过与与非门和或门的处理，再使用其输出经过与门的处理就完成了异或的操作：这种多层感知机就是一种二层神经网络，层数的算法和感知机的层数算法一致，只是表示法不同：先解释结构，x1 x2的位置称为输入层，y所在位置称为输出层，两层之间的不管有多少层都是中间层或者隐藏层，不同点在于门的运算过程都跑到边上去了，回顾感知机的实现，所谓门的实现就是不同权重的调整，所以边上既包含了参数也包含了固定的传播算法，节点或者称神经元上包含了输出输入的结果，但也包含了一些其他处理过程，如激活函数一般而言，“朴素感知机”是指单层网络，指的是激活函数使用了阶跃函数(后述）的模型。“多层感知机”是指神经网络，即使用sigmoid函数(后述）等平滑的激活函数的多层网络。传播算法和激活函数首先，回顾下感知机中的权重和偏置和感知机的实现：简单来说感知机的实现就是该公式的实现，其中b是被称为偏置的参数，用于控制神经元被激活的容易程度；而w1和w2是表示各个信号的权重的参数，用于控制各个信号的重要性。 而b+w1x1+w2x2就是层和层之间的传播算法，从下图可以看出，偏置也可作为神经元为1而权重为b进行传播：解决了传播问题，再看一眼公式可以发现0和1的判断是依赖于b+w1x1+w2x2的输出的，而这0和1的判断规则也不一定是以公式中的方式出现的，是可以记作函数h(x)来表示的，这个函数就是激活函数，激活函数的输出才是神经元最后表示的输出。简而言之：激活函数，是用来决定何时激活信号的函数阶跃函数就从感知机使用的激活函数来看：该函数称为阶跃函数；实现和图像如下，注意到这并不是一个平滑的函数def step_function(x): y = x > 0 return y.astype(int) x = np.array([-1.0,1.0,2.0]) print(x) print(x>0) y = step_function(x) print(type(y)) #[-1. 1. 2.] #[False True True] #<class 'numpy.ndarray'> import matplotlib.pyplot as plt x = np.arange(-5.,5.,0.1) y = step_function(x) plt.plot(x,y) plt.ylim(-0.1,1.1) plt.show()Sigmoid函数Sigmoid函数是神经网络经常使用一个函数，公式，实现和图像如下：def sigmoid(x): return 1/(1+np.exp(-x)) x = np.arange(-5.,5.,0.1) y = sigmoid(x) plt.plot(x,y) plt.ylim(-0.1,1.1) plt.show()函数对比上面说到，阶跃函数并不是平滑的函数，而sigmoid是，然后二者都是非线性函数，这里先给出结论：神经网络的激活函数必须使用非线性函数。原因是由线性函数的性质造成的，如果把中间层合并在一起，所做的处理更加类似于一个超级复杂的复合函数，但是线性函数再怎么复合，都可以找到一个参数来等价的表示这些复合的操作，对于中间层来说就是可以用一层来代表N层，这样中间层的存在意义就没有了，神经网络也就无法发挥多层网络带来的优势。下图为简单的举例：ReLU函数sigmoid函数很早就开始被使用了，而最近则主要使用ReLU（Rectified Linear Unit）函数。公式，实现和图像如下：def relu(x): return np.maximum(0,x) x = np.arange(-5.,5.,0.1) y = relu(x) plt.plot(x,y)神经网络和矩阵矩阵乘法回顾一下numpy的相关操作，shape表示数组形状，dim表示数组维数；多维数组可以代表矩阵，横为行竖为列；然后简单了解下矩阵乘法：代码如下：A=np.array([[1,2,3],[4,5,6]]) B=np.array([[1,2],[3,4],[5,6]]) print(A.shape,B.shape) C=np.dot(A,B) print(C) print(C.shape) #(2, 3) (3, 2) #[[22 28] # [49 64]] #(2, 2)在矩阵的乘积运算中，对应维度的元素个数要保持一致 ：神经网络和矩阵从下图可见，之前所说层与层之间的传播算法，可以统统用矩阵乘法表示。而这一层接一层的传播就是所谓的（前）正向传播。前向传播前向传播用语言表达比较繁琐，直接看图和代码则比较方便，无非是一层接一层的矩阵计算和激活函数计算：输出层的激活函数激活函数已经了解，但是输出层的激活函数稍微有些不同，但也只是输出函数类型的问题。首先 ，对于机器学习解决的问题有两种，一种是所谓分类问题，就好比数字识别，从图片判断是0-9的哪一个，还有一种是回归问题，好比从历史数据判断今天地铁的客流量，一个是选择题，一个是填空题。在输出层，对于选择题使用softmax函数，对于回归问题使用恒等函数。顺带一提，输出层神经元的数量在分类问题中取决于分类的个数。恒等函数恒等函数图像如下，其实就是将输入原样输出，这里便不列出代码identity_function了SoftMax函数对于分类使用的SoftMax函数则有一些说法，先看公式和基本实现:首先这涉及了前一层的所有输入，类比sql语句，这便成为全连接，意思是每一个输出都完整的接受了上一个输入层中的所有数据；而计算时用所有输入做分母，其中一个做分子，这又有点概率的意思，事实上也是这样最后，SoftMax函数一般只用于学习阶段，因为在推理过程中，一般而言只把输出值最大的神经元所对应的类别作为识别结果，使用softmax会造成额外的计算量下面是代码实现：def softmax(a): exp_a = np.exp(a) sum_exp_a=np.sum(exp_a) y=exp_a/sum_exp_a return y a = np.array([0.3,2.9,4.0]) print(softmax(a)) #下两行会引发溢出异常 a = np.array([1010,1000,990]) print(softmax(a)) 这样是用问题的，因为对于指数函数有溢出的风险，但如果按如下步骤进行改造则可以解决这个问题：其中 C'一般取a_n这个数组中的最大值取反，这样削减一下大小，就不会有溢出的问题了：# softmax的改进变形 def softmax(a): c = np.max(a) exp_a = np.exp(a-c) #溢出对策 sum_exp_a = np.sum(exp_a) y = exp_a /sum_exp_a return y a = np.array([1010,1000,990]) print(softmax(a)) #[9.99954600e-01 4.53978686e-05 2.06106005e-09]总结和阶段性实现目前，一个支持前向传播的神经网络便可以实现，但是目前只能假设一个权重参数（因为没有学习阶段的支持），然后通过矩阵乘法和中间层和输出层的激活函数来实现前向传播。所以目前以MNIST数据集为例来实现一个简易的前向传播测试数据和训练数据机器学习中，一般将数据分为训练数据（监督数据 ）和测试数据两部分来进行学习和实验等。首先，使用训练数据进行学习，寻找最优的参数；然后，使用测试数据评价训练得到的模型的实际能力在处理测试数据过程中，如果发现效果很好，则说明该模型的泛化能力很好，即面对未知数据表现很好；反面的，如果只对训练数据表现好，则说明泛化能力差，这种现象称为过拟合MNIST数据集的读取读取该数据集需要一些工具函数的辅助，这里直接使用《深度学习入门》中提供的工具实现：# 这里的路径按本地实际来 from dfs.dataset.mnist import load_mnist import sys,os sys.path.append(os.pardir) (x_train, t_train), (x_test, t_test) = load_mnist(flatten=True, normalize=False,one_hot_label=False) # 输出各个数据的形状 print(x_train.shape) # (60000, 784) # 输出各个数据的形状 print(t_train.shape) # (60000,) print(x_test.shape) # (10000, 784) print(t_test.shape) # (10000,)train代表训练数据，test表示测试数据x代表输入的图片，t代表最终的标签（在one_hot_label是false时会在显示标签0-9,在true时会返回一个10位数组，在标签数字对应的位置是1，其余为0）from PIL import Image def img_show(img): pil_img = Image.fromarray(np.uint8(img)) pil_img.show() img = x_train[0] label = t_train[0] print(label) # 5 print(img.shape) # (784,) img = img.reshape(28, 28) # 把图像的形状变成原来的尺寸 print(img.shape) # (28, 28) img_show(img)该数据集中，输入的是一堆图片和图片对应的标签，在这里将图片展开为一维数组(flatten=True)，并使用原始像素信息，而不是将0-255的颜色值压缩到0-1的范围(normalize=false)，而将数据压缩到0-1范围这种操作称为标准化或正交化one-hot表示是仅正确解标签为1，其余皆为0的数组，就像[0,0,1,0,0,0,0,0,0,0]这样。当one_hot_label为False时，只是像 7、 2这样简单保存正确解标签；当 one_hot_label为 True时，标签则保存为one-hot表示。实现首先，由于并没有学习这一步骤，所以会给所有权重和偏置设置符合正态分布的默认值；然后因为传播的实现就是矩阵乘法的实现，现在假设使用的是3层神经网络，则三个权重矩阵的形状要符合矩阵乘法的要求，所以假设矩阵形状如下：最后，由于一个图片一个图片处理过于缓慢，而矩阵乘法却能吃到SIMD的优化福利，所以会使用批处理的方法，即如下图所示，一次性打包处理100张图像，为此，可以把x的形状改为100 × 784，将100张图像打包作为输入数据：如此实现如下，首先是神经网络的实现和使用import numpy as np def softmax(a): c = np.max(a) exp_a = np.exp(a-c) #溢出对策 sum_exp_a = np.sum(exp_a) y = exp_a /sum_exp_a return y def sigmoid(x): return 1/(1+np.exp(-x)) class SigMoidLayer: def __init__(self,w,b,s="SigMoid"): self.layer_name = s self.w=w self.b = b def forward(self,x): a=np.dot(x,self.w)+self.b return sigmoid(a) class OutputLayer2: def __init__(self,w,b,s="Ouput"): self.layer_name = s self.w=w self.b = b def forward(self,x): a=np.dot(x,self.w)+self.b return softmax(a) class TestMNIST: # 正常情况下，权重或者偏置都应该是学习的结果，这里直接赋值只是一种假设 def __init__(self): self.network_w = {} self.network_b = {} self.network_w['L1'] = np.random.randn(784,50) self.network_b['L1'] = 1 self.network_w['L2'] = np.random.randn(50,100) self.network_b['L2'] = 1 # 输出层 self.network_w['L3'] = np.random.randn(100,10) self.network_b['L3'] = 1 self.layers=['L1','L2'] def forward(self,x): tmp=x for l in self.layers: tmp = SigMoidLayer(self.network_w[l],self.network_b[l],l).forward(tmp) print(l,tmp) tmp = OutputLayer2(self.network_w['L3'],self.network_b['L3'],'L3').forward(tmp) return tmp nett = TestMNIST() #测试代码 #nett.forward(np.random.randn(100,784)) batch_size=100 for i in range(0,len(x_test),batch_size): x_batch = x_test[i:i+batch_size] #print(x_batch.shape) y_batch = nett.forward(x_batch) p=np.argmax(y_batch,axis=1) print(p) #[8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 # 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 # 8 8 8 8 8 8 6 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8] 可见结果并不好，都只会输出8，即使如此也要测一下准确率：accuracy_cnt=0 for i in range(0,len(x_test),batch_size): x_batch = x_test[i:i+batch_size] y_batch = nett.forward(x_batch) p=np.argmax(y_batch,axis=1) accuracy_cnt += np.sum(p == t_test[i:i+batch_size]) print("Accuracy:" + str(float(accuracy_cnt) / len(x_test))) #Accuracy:0.0965太好了，基本都对不上。
通于一而万事毕，无心得而鬼神服。本文以数组加法为例，简述CUDA的编程模型基本概念确保了解《CUDA的线程结构》中的内容再进行下去，因为对于线程结构本章会省略一部分描述编程模型首先，编程模型是写程序时可被控制的部分，通过编程模型可以控制计算设备的行为；其次，编程模型对于编译器和库函数来说是一层通信抽象；再次，对于GPU来说，编程模型可分为核函数内存管理线程管理流最后，对于Cuda编程，可抽象出几个过程领域层（业务主导要实现的目标）逻辑层（基于逻辑组织线程代码）硬件层（思考如何将线程映射到机器上，以及思考如何提升性能）接下来，会从以下方面说明Cuda的编程结构内存线程核函数启动编写验证错误处理内存管理首先应该理解异构环境下内存的样貌：通过PCIe总线分割和通信的两种内存：主机内存和设备内存一个异构环境，通常有多个CPU多个GPU，他们都通过PCIe总线相互通信，也是通过PCIe总线分隔开的。所以我们要区分一下两种设备的内存：主机：CPU及其内存设备：GPU及其内存注意这两个内存从硬件到软件都是隔离的（CUDA6.0 以后支持统一寻址），我们目前先不研究统一寻址，我们现在还是用内存来回拷贝的方法来编写调试程序，以巩固大家对两个内存隔离这个事实的理解。CUDA的内存管理是C风格的，列表如下：标准C函数CUDA C 函数说明malloccudaMalloc内存分配memcpycudaMemcpy内存复制memsetcudaMemset内存设置freecudaFree释放内存另外，CUDA中设备内存也是分层次的，分为单网格中公用的全局内存和单线程块中公用的共享内存，大致如图：线程管理在《CUDA的线程结构》中，曾经用草图展示过cuda线程的组织形式，这里用更好的图展示一下：另外，不同线程块的线程物理隔离无法相互影响，同一个线程块中的线程可以实现以下两种操作：同步共享内存最后再重申一下《CUDA的线程结构》中提过的，如何标识一个线程：blockIdx（线程块在线程网格内的位置索引）threadIdx（线程在线程块内的位置索引）这两个内置结构体基于 uint3 定义，包含三个无符号整数xyz的结构还有blockIdx中三个字段的范围和threadIdx中三个字段的范围：gridDimblockDimdim3类型(基于uint3定义的数据结构)的变量，也包含三个字段xyz另外，线程块在3个轴上的个数和一个线程块的线程数是有限制的：线程块在xyz三轴上分别最多：2^31-1,65535,65535一个线程块线程总数最多1024，在三维上最多:1024,1024,64核函数之前hello world时已经用到了核函数，但是没有给一个明确定义，这里给出：核函数就是在CUDA模型上诸多线程中运行的那段串行代码，这段代码在设备上运行，用NVCC编译，产生的机器码是GPU的机器码，所以我们写CUDA程序就是写核函数第一步我们要确保核函数能正确的运行产生正切的结果第二步优化CUDA程序的部分，无论是优化算法，还是调整内存结构，线程结构都是要调整核函数内的代码，来完成这些优化的。接下来将以数组加法为例来实现和使用核函数，以及应用一些内存管理的函数单机版数组加法无需多言，直接上代码，明确我们要干什么：/** * 单机版:两数组相加 * */ #include<iostream> const double a=1.23; const double b=3.45; const double c=4.68; const int N = 100000000; const double EPSILON = 1.0e-15; // 函数1: 数组相加(数组默认是动态分配的) void add(const double *x,const double *y,double *z,const int N) { for (int i=0;i<N;i++) { z[i] = x[i]+y[i]; } } // 函数2: 检查结果 void check(const double *z,const int N) { bool has_error = false; for (int i=0;i<N;i++) { if (fabs(z[i]-c) > EPSILON) { has_error = true; break; } } printf("%s\n",has_error?"Has errors":"No errors"); } int main() { //创建数组xyz double* x = new double[N]; double* y = new double[N]; double* z = new double[N]; //初始化数组 for (int i=0;i<N;i++) { x[i] = a; y[i] = b; z[i] = 0; } //计算数组 add(x,y,z,N); //检查结果 check(z,N); //删除数组指针 delete[] x; delete[] y; delete[] z; return 0 ; }CUDA程序和内存管理先看代码，比起hello world，这是一个完整的CUDA程序，它包含了下图所示的结构：/** * cuda版:两数组相加 * */ #include<iostream> const double a=1.23; const double b=3.45; const double c=4.68; const int N = 100000000; const int M = sizeof(double)*N; const double EPSILON = 1.0e-15; // 函数1: 数组相加(数组默认是动态分配的) __global__ void add(const double *x,const double *y,double *z) { const int n = blockDim.x*blockIdx.x + threadIdx.x; z[n] = x[n] + y[n]; } // 函数2: 检查结果 void check(const double *z,const int N) { bool has_error = false; for (int i=0;i<N;i++) { if (fabs(z[i]-c) > EPSILON) { has_error = true; break; } } printf("%s\n",has_error?"Has errors":"No errors"); } int main() { //主机: 创建数组xyz double* x = new double[N]; double* y = new double[N]; double* z = new double[N]; //初始化数组 for (int i=0;i<N;i++) { x[i] = a; y[i] = b; z[i] = 0; } //主机: 创建设备数组xd,yd double* xd; double* yd; cudaMalloc((void **)&xd,M); cudaMalloc((void **)&yd,M); //主机: 复制主机数据到设备 cudaMemcpy(xd,x,M,cudaMemcpyHostToDevice); cudaMemcpy(yd,y,M,cudaMemcpyHostToDevice); //主机: 指定设备运行参数 //const dim3 block_size(128,1,1); const int block_size=128; const dim3 grid_size(N/block_size,1,1); //设备: 计算数组 add<<<grid_size,block_size>>>(xd,yd,yd); //主机: 复制设备结果到主机 cudaMemcpy(z,yd,M,cudaMemcpyDeviceToHost); //检查结果 check(z,N); //删除数组指针 cudaFree(xd); cudaFree(yd); delete[] x; delete[] y; delete[] z; return 0 ; }先把核函数的部分放在一边，代码里有以下几点可以注意：隐式设备初始化：了解过opengl或d3d程序的可以明显感觉，CUDA是没有显式初始化设备(即GPU)的函数。在第一次调用一个和设备管理及版本查询功能无关的运行时 API函数时，设备将自动初始化。内存管理一节已经大致了解了内存管理的相关函数，这里就是应用的地方，需要注意的是cudaMalloc和cudaMemcpycudaMalloc的第一个参数是指针的指针cudaMemcpy的最后一个参数规定了数据传输的方式（从名称就可以看出）CUDA程序和核函数然后来看核函数，可分为三个方面：启动核函数编写核函数验证核函数首先，CPU是一个控制者，运行核函数，要从CPU发起，第一，核函数是异步的，第二，如《CUDA的线程结构》所说，三个尖括号<<<grid,block>>>内是对设备代码执行的线程结构的配置：网格的布局（线程块）线程块的布局（线程）其次，同步函数有两种cudaDeviceSynchronize显式同步涉及到设备端不执行完，主机没办法进行，比如内存拷贝函数，隐式同步然后看到核函数的写法，cuda中不仅有核函数，还有一些其他类型，通过不同的标识表明，如核函数是__global__：限定符执行调用备注global设备端执行可以从主机调用也可以从计算能力3以上的设备调用必须有一个void的返回类型device设备端执行设备端调用 host主机端执行主机调用可以省略有些函数可以同时定义为 device 和 host ，这种函数可以同时被设备和主机端的代码调用，主机端代码调用函数很正常，设备端调用函数与C语言一致，但是要声明成设备端代码，告诉nvcc编译成设备机器码，同时声明主机端设备端函数，那么就要告诉编译器，生成两份不同设备的机器码。其他方面核函数基本和普通函数一致，但有一些限制：除了统一寻址，只能访问设备内存必须有void返回类型不支持可变数量的参数不支持静态变量显示异步行为不可作为成员，但可以包装在动态并行的前提下，可以递归调用另外设备函数可用__noinline__ __forceinline__来作为内联函数和非内联函数CUDA小技巧，当我们进行调试的时候可以把核函数配置成单线程的：下面是一些例子：函数的使用斐波那契数列，本程序意义不大，旨在展示__device__ __host__的共用斐波那契数列，本程序意义不大，旨在展示递归的使用/** * 函数的使用 */ #include<iostream> const double a=1.23; const double b=3.45; const double c=4.68; const int N = 100000000; const int M = sizeof(double)*N; const double EPSILON = 1.0e-15; // 函数1: 数组相加(数组默认是动态分配的) //TODO: 注意到三种写法效率是一样的 __device__ double addfunc1(double x,double y) { return x+y; } __device__ double addfunc2(double x,double y,double *z) { *z = x+y; } __device__ double addfunc3(double x,double y,double &z) { z= x+y; } __global__ void add(const double *x,const double *y,double *z) { const int n = blockDim.x*blockIdx.x + threadIdx.x; if (n<N) { z[n] = x[n] + y[n]; //z[n] = addfunc1(x[n],y[n]); //addfunc2(x[n],z[n],&z[n]); //addfunc3(x[n],y[n],z[n]); } }/** * 斐波那契数列，本程序意义不大，旨在展示__device__ __host__的共用 * */ #include<iostream> const int N = 10; const int M = sizeof(int)*N; // 函数1: 斐波那契数列 __device__ __host__ void fibN(int* A, const int N) { if (N<0) return; if (A[N]==-1) { if (N>1) { fibN(A,N-1); fibN(A,N-2); A[N] = A[N-1] + A[N-2]; } if (N<=1) { A[N]=N; } } } __global__ void fibNN(int* A, const int N) { fibN(A, N); } // 函数2: 检查结果 //0、1、1、2、3、5、8、13、21、34 void check(const int *A,const int N) { for (int i=0;i<N;i++) { std::cout<<A[i]<<" "; } std::cout<<std::endl; } int main() { //主机: 创建数组xyz int* x = new int[N]; //初始化数组 for (int i=0;i<N;i++) { x[i] = -1; } //主机: 创建设备数组xd,yd int* xd; cudaMalloc((void **)&xd,M); //主机: 复制主机数据到设备 cudaMemcpy(xd,x,M,cudaMemcpyHostToDevice); //主机: 指定设备运行参数 const dim3 block_size(128,1,1); const dim3 grid_size(128,1,1); //设备: 计算数组 fibNN<<<grid_size,block_size>>>(xd,N-1); //主机: 复制设备结果到主机 cudaMemcpy(x,xd,M,cudaMemcpyDeviceToHost); //检查结果 check(x,N); for (int i=0;i<N;i++) { x[i] = -1; } fibN(x,N-1); check(x,N); //删除数组指针 cudaFree(xd); delete[] x; return 0 ; }/** * cuda版:两数组相加 * */ #include<iostream> const double a=1.23; const double b=3.45; const double c=4.68; const int N = 100000000; const int M = sizeof(double)*N; const double EPSILON = 1.0e-15; // 函数1: 数组相加(数组默认是动态分配的) __global__ void add(const double *x,const double *y,double *z) { const int n = blockDim.x*blockIdx.x + threadIdx.x; if (n<N) { z[n] = x[n] + y[n]; } } // 函数2: 检查结果 void check(const double *z,const int N) { bool has_error = false; for (int i=0;i<N;i++) { if (fabs(z[i]-c) > EPSILON) { has_error = true; break; } } printf("%s\n",has_error?"Has errors":"No errors"); } int main() { //主机: 创建数组xyz double* x = new double[N]; double* y = new double[N]; double* z = new double[N]; //初始化数组 for (int i=0;i<N;i++) { x[i] = a; y[i] = b; z[i] = 0; } //主机: 创建设备数组xd,yd double* xd; double* yd; cudaMalloc((void **)&xd,M); cudaMalloc((void **)&yd,M); //主机: 复制主机数据到设备 cudaMemcpy(xd,x,M,cudaMemcpyHostToDevice); cudaMemcpy(yd,y,M,cudaMemcpyHostToDevice); //主机: 指定设备运行参数 //const dim3 block_size(128,1,1); const int block_size=128; const dim3 grid_size(N/block_size,1,1); std::cout << "Size: " << N/block_size <<"=="<<block_size <<std::endl; //设备: 计算数组 //TODO: 注意到这里使用了会超过N的执行参数 add<<<800000,block_size>>>(xd,yd,yd); //主机: 复制设备结果到主机 cudaMemcpy(z,yd,M,cudaMemcpyDeviceToHost); //检查结果 check(z,N); //删除数组指针 cudaFree(xd); cudaFree(yd); delete[] x; delete[] y; delete[] z; return 0 ; }/** * 斐波那契数列，本程序意义不大，旨在展示递归的使用 * */ #include<iostream> const int N = 10; const int M = sizeof(int)*N; // 函数1: 斐波那契数列 __global__ void fibN(int* A, const int N) { if (N<0) return; if (A[N]==-1) { if (N>1) { fibN<<<1,1>>>(A,N-1); fibN<<<1,1>>>(A,N-2); A[N] = A[N-1] + A[N-2]; } if (N<=1) { A[N]=N; } } } // 函数2: 检查结果 //0、1、1、2、3、5、8、13、21、34 void check(const int *A,const int N) { for (int i=0;i<N;i++) { std::cout<<A[i]<<" "; } std::cout<<std::endl; } int main() { //主机: 创建数组xyz int* x = new int[N]; //初始化数组 for (int i=0;i<N;i++) { x[i] = -1; } //主机: 创建设备数组xd,yd int* xd; cudaMalloc((void **)&xd,M); //主机: 复制主机数据到设备 cudaMemcpy(xd,x,M,cudaMemcpyHostToDevice); //主机: 指定设备运行参数 const dim3 block_size(128,1,1); const dim3 grid_size(128,1,1); //设备: 计算数组 fibN<<<grid_size,block_size>>>(xd,N-1); //主机: 复制设备结果到主机 cudaMemcpy(x,xd,M,cudaMemcpyDeviceToHost); //检查结果 check(x,N); //删除数组指针 cudaFree(xd); delete[] x; return 0 ; } //0 1 1 0 -2 -2 -2 -4 -2 -6 结果错误 但是可见可以这么用验证和错误检测所有编程都需要对错误进行处理，早期的编码错误，编译器会帮我们搞定，内存错误也能观察出来，但是有些逻辑错误很难发现，甚至到了上线运行时才会被发现，而且有些厉害的bug复现会很难，不总出现，但是很致命，而且CUDA基本都是异步执行的，当错误出现的时候，不一定是哪一条指令触发的，这一点非常头疼；这时候我们就需要对错误进行防御性处理了错误检测宏如下：#pragma once #include <cstdio> #define CHECK(call) \ do \ { \ const cudaError_t error_code = call; \ if (error_code != cudaSuccess) \ { \ printf("CUDA Error:\n"); \ printf(" File: %s\n", __FILE__); \ printf(" Line: %d\n", __LINE__); \ printf(" Error code: %d\n", error_code); \ printf(" Error text: %s\n", \ cudaGetErrorString(error_code)); \ exit(1); \ } \ } while (0) cudaGetErrorString可把error_code转为相应的提示信息用这个宏函数包装大部分CUDA运行时API函数。有一个例外是cudaEventQuery()函数，因为它很有可能返回cudaErrorNotReady,但又不代表程序出错了#include "error.cuh" #include <cmath> #include <cstdio> const double EPSILON = 1.0e-15; const double a = 1.23; const double b = 2.34; const double c = 3.57; void __global__ add(const double *x, const double *y, double *z, const int N); void check(const double *z, const int N); int main(void) { const int N = 100000000; const int M = sizeof(double) * N; double *h_x = (double*) malloc(M); double *h_y = (double*) malloc(M); double *h_z = (double*) malloc(M); for (int n = 0; n < N; ++n) { h_x[n] = a; h_y[n] = b; } double *d_x, *d_y, *d_z; CHECK(cudaMalloc((void **)&d_x, M)); CHECK(cudaMalloc((void **)&d_y, M)); CHECK(cudaMalloc((void **)&d_z, M)); //TODO:此处错误复制方向搞反了 CHECK(cudaMemcpy(d_x, h_x, M, cudaMemcpyDeviceToHost)); CHECK(cudaMemcpy(d_y, h_y, M, cudaMemcpyDeviceToHost)); const int block_size = 128; const int grid_size = (N + block_size - 1) / block_size; add<<<grid_size, block_size>>>(d_x, d_y, d_z, N); CHECK(cudaMemcpy(h_z, d_z, M, cudaMemcpyDeviceToHost)); check(h_z, N); free(h_x); free(h_y); free(h_z); CHECK(cudaFree(d_x)); CHECK(cudaFree(d_y)); CHECK(cudaFree(d_z)); return 0; } void __global__ add(const double *x, const double *y, double *z, const int N) { const int n = blockDim.x * blockIdx.x + threadIdx.x; if (n < N) { z[n] = x[n] + y[n]; } } void check(const double *z, const int N) { bool has_error = false; for (int n = 0; n < N; ++n) { if (fabs(z[n] - c) > EPSILON) { has_error = true; } } printf("%s\n", has_error ? "Has errors" : "No errors"); }用上述方法不能捕捉调用核函数的相关错误，因为核函数不返回任何值(回顾一下，核函数必须用void修饰)。有一个方法可以捕捉调用核函数可能发生的错误#include "error.cuh" #include <cmath> #include <cstdio> const double EPSILON = 1.0e-15; const double a = 1.23; const double b = 2.34; const double c = 3.57; void __global__ add(const double *x, const double *y, double *z, const int N); void check(const double *z, const int N); int main(void) { const int N = 100000000; const int M = sizeof(double) * N; double *h_x = (double*) malloc(M); double *h_y = (double*) malloc(M); double *h_z = (double*) malloc(M); for (int n = 0; n < N; ++n) { h_x[n] = a; h_y[n] = b; } double *d_x, *d_y, *d_z; CHECK(cudaMalloc((void **)&d_x, M)); CHECK(cudaMalloc((void **)&d_y, M)); CHECK(cudaMalloc((void **)&d_z, M)); CHECK(cudaMemcpy(d_x, h_x, M, cudaMemcpyHostToDevice)); CHECK(cudaMemcpy(d_y, h_y, M, cudaMemcpyHostToDevice)); //TODO:错误在线程块的大小太大 const int block_size = 1280; const int grid_size = (N + block_size - 1) / block_size; add<<<grid_size, block_size>>>(d_x, d_y, d_z, N); CHECK(cudaGetLastError()); CHECK(cudaDeviceSynchronize()); CHECK(cudaMemcpy(h_z, d_z, M, cudaMemcpyDeviceToHost)); check(h_z, N); dim3 a; free(h_x); free(h_y); free(h_z); CHECK(cudaFree(d_x)); CHECK(cudaFree(d_y)); CHECK(cudaFree(d_z)); return 0; } void __global__ add(const double *x, const double *y, double *z, const int N) { const int n = blockDim.x * blockIdx.x + threadIdx.x; if (n < N) { z[n] = x[n] + y[n]; } } void check(const double *z, const int N) { bool has_error = false; for (int n = 0; n < N; ++n) { if (fabs(z[n] - c) > EPSILON) { has_error = true; } } printf("%s\n", has_error ? "Has errors" : "No errors"); }CHECK(cudaGetLastError()); CHECK(cudaDeviceSynchronize());其中，第一条语句的作用是捕捉第二个语句之前的最后一个错误，第二条语句的作用是同步主机与设备。之所以要同步主机与设备，是因为核函数的调用是异步的，即主机发出调用核函数的命令后会立即执行后面的语句，不会等待核函数执行完毕。
万物一齐，而百兽率舞。本文简略介绍CUDA编程的相关背景异构计算首先，异构计算就是CPU+X，这个X可以是GPU，也可以是NPU，FPGA，DSP等等，在CPU上运行的代码称为Host主机代码，X上运行的代码称为Device设备代码其次，异构计算和同构计算的一个不同点在于，异构的处理要比同构复杂，因为涉及同构下自动处理的计算，控制，传输都要人工进行干预GPUCPU vs GPU先看cpu的架构和gpu架构的不同：先看几个共同的部分：控制器Control缓存器Cache随机存取存储器DRAM总线PCIe Bus算术逻辑单元ALU不需要过多的计组原理知识，直观来看，GPU弱化了控制器的能力，大大堆砌了ALU的数量，且为一组ALU（SM，简单理解为线程束）分配了单独的Control和cache，所以直觉上，对了逻辑简单，数据量大的任务，GPU更高效低并行逻辑复杂的程序适合用CPU高并行逻辑简单的大数据计算适合GPUCPU和GPU线程的区别：CPU线程是重量级实体，操作系统交替执行线程，线程上下文切换花销很大GPU线程是轻量级的，GPU应用一般包含成千上万的线程，多数在排队状态，线程之间切换基本没有开销。CPU的核被设计用来尽可能减少一个或两个线程运行时间的延迟，而GPU核则是大量线程，最大幅度提高吞吐量GPU简介对于Nvida，GPU分为4个系列：TegraGeforceQuadroTesla针对不同的应用场景，比如Tegra用于嵌入式，Geforce是平时打游戏用到，Tesla主要用于计算，但目前看Geforce和Tesla貌似也没那么分得很清楚，只是某些系列上有容错的硬件部分GPU也存在版本号，即其计算能力，但这和性能不是成正比的，只代表架构的更新，对于CUDA来说，新版的cuda只支持最近几个版本的GPU，参考：https://en.wikipedia.org/wiki/CUDA#GPUs_supported真正评价GPU能力的分为两种，一种是计算性能指标，一种是容量指标：计算性能主要看FLOPS，一般单位在T级别（TFLOPS），称每秒浮点数运算次数，分单精度双精度，根据不同的架构，双精度一般是单精度的1/N容量指标上有：显存，带宽，核心数量等，特别的如果一个CUDA程序爆显存了，则不使用统一内存的情况下是无法运行程序的CUDACUDA平台不是单单指软件或者硬件，而是建立在Nvidia GPU上的一整套平台，并扩展出多语言支持对于CUDA API来说，分为driver和运行时两种API，二者互斥不能同时使用，性能没啥差距driver API 低级，复杂，灵活，一般用于配合其他语言调用CUDA运行时简单，基于driver API，稍微高级一点CUDA编程环境架构如下：nvccnvcc之于cuda相当于普通c++之于g++。因为cuda程序中host和device的代码是混在一起的，nvcc会分离这两个部分分别处理：其中，PTX相当于CUDA的汇编语言nvdia-smi可以通过nvidia-smi(Nvidia's system management interface)程序检查与设置设备。它包含在CUDA开发工具套装内。CUDA的程序结构一般CUDA程序分成下面这些步骤：分配GPU内存拷贝内存到设备调用CUDA内核函数来执行计算把计算完成数据拷贝回主机端内存销毁参考https://face2ai.com/program-blog/#GPU%E7%BC%96%E7%A8%8B%EF%BC%88CUDA%EF%BC%89
万物化作，萌区有状，盛衰之杀，变化之流也。本文以helloworld为例，简述CUDA的线程结构基本概念即使是HelloWorld，Cuda程序也引入了一些不同的概念，确保了解《简略CUDA和GPU相关概念》中的内容在进行下去从HelloWorld开始下面是一个最简单的HelloWorld程序：/** 运行命令: cl ./hello_host_only.cpp **/ #include <cstdio> using namespace std; int main() { printf("Host: Hello wolrd!"); return 0; }在windows可在x64 Native Tools Command Prompt for VS 2022用cl来编译它只有Host代码的CUDA版本HelloWolrd/** 运行命令: nvcc -arch=sm_62 hello_host_only.cu sm_62的来源: https://docs.nvidia.com/cuda/archive/12.9.0/cuda-compiler-driver-nvcc/index.html https://en.wikipedia.org/wiki/CUDA#GPUs_supported */ #include <cstdio> using namespace std; int main() { printf("Host: Hello wolrd!"); return 0; }代码完全没有变化，只有文件后缀名发生了变化但这里性质发生了变化，代码分为了Host代码部分和Device代码部分由nvcc分别处理nvcc指定运算架构注意到注释部分，nvcc可以指定虚拟和真实架构的编号，可以通过链接查看本地设备可以使用的标识如果需要匹配多个架构，可以使用多个-gencode -arch 虚拟架构 -code 真实架构来编译产生不同的输出，最后真实架构不能落后与虚拟架构使用CUDA线程从这里正式开始添加Device代码：/** 运行命令: nvcc XX.cu */ #include<cuda.h> #include<cuda_runtime.h> #include<cstdio> __global__ void hello_from_gpu() { printf("Device: Hello World."); } int main() { const int grid_size=1; const int block_size=1; hello_from_gpu<<<grid_size,block_size>>>(); cudaDeviceSynchronize(); return 0; }注意四点：include包含了两个cuda相关头文件，这不是必须要添加的，nvcc会自动引入，其中包含了cstdlib的内容，但没有包含cstdio的内容，最后device程序是不认iostream的__global___就是device代码，即核函数的标志cudaDeviceSynchronize的功能是同步host和device，也就是说gpu上的代码和cpu上的代码是异步运行的<<<grid_size,block_size>>>指定了用多少线程运行该核函数使用多线程的CUDA程序再进一步改进程序/** 运行命令: nvcc XX.cu */ #include<cuda.h> #include<cuda_runtime.h> #include<cstdio> __global__ void hello_from_gpu() { const int bid = blockIdx.x; const int tid = threadIdx.x; printf("Device: Hello World from block %d and thread %d\n",bid,tid); } int main() { const int grid_size=2; const int block_size=4; hello_from_gpu<<<grid_size,block_size>>>(); cudaDeviceSynchronize(); return 0; }CUDA将线程组织为网格，网格被划分为线程块，即网格由线程块组成，线程块由线程组成，三者之间并不是简单的包含关系（即一般论的包含关系，网格有N个线程块，每个线程块又有M个线程），但是目前示例中确实只是单纯的包含关系，可简单认为有网格两个线程块，每个线程块四个线程共八个线程下图展示了多维线程的布局，这也是网格的真正分布：可以认为网格是一个三维立方体，被划分出的小格就是线程块，而小格又可以按三维继续分割，分割出的就是线程。像上面简单的以2和4指定的线程块个数和线程个数，三维上可以分别认为是(2,1,1)和(4,1,1)，只按x轴延伸的2个大方块，每个大方块又是由按x轴延展的4个小方块组成多维线程的CUDA程序既然有了三维建构，每个线程块和线程自然有其三维坐标，按照简单的数学逻辑也可以转换为1维坐标，这里直接把这些概念在代码中展示：/** 运行命令: nvcc XX.cu */ #include<cuda.h> #include<cuda_runtime.h> #include<cstdio> __global__ void hello_from_gpu() { const int bid = blockIdx.x; const int tid_x = threadIdx.x; const int tid_y = threadIdx.y; const int tid_z = threadIdx.z; //1D坐标: const int bid_1D = blockIdx.x+blockIdx.y*gridDim.x+blockIdx.z*gridDim.x*gridDim.y; const int tid_1D = threadIdx.x + blockDim.x*threadIdx.y + blockDim.x*blockDim.y*threadIdx.z; printf("Device(%d)(%d): Hello World from block %d and thread (%d-%d-%d)\n",bid_1D,tid_1D,bid,tid_x,tid_y,tid_z); } int main() { const int grid_size=2; const dim3 block_size(3,3,3); hello_from_gpu<<<grid_size,block_size>>>(); cudaDeviceSynchronize(); return 0; }注意到：dim3类似于vec3，表3维向量gridDim和blockDim是dim3固有变量，内容是当初规定的各个轴的上限blockIdx和threadIdx也是dim3固有变量，包含了当前线程的坐标另外，线程块在3个轴上的个数和一个线程块的线程数是有限制的：线程块在xyz三轴上分别最多：2^31-1,65535,65535一个线程块线程总数最多1024，在三维上最多:1024,1024,64下面是超过了限制的程序：/** 运行命令: nvcc XX.cu */ #include<cuda.h> #include<cuda_runtime.h> #include<cstdio> __global__ void hello_from_gpu() { const int bid = blockIdx.x; const int tid_x = threadIdx.x; const int tid_y = threadIdx.y; const int tid_z = threadIdx.z; //1D坐标: const int bid_1D = blockIdx.x+blockIdx.y*gridDim.x+blockIdx.z*gridDim.x*gridDim.y; const int tid_1D = threadIdx.x + blockDim.x*threadIdx.y + blockDim.x*blockDim.y*threadIdx.z; printf("Device(%d)(%d): Hello World from block %d and thread (%d-%d-%d)\n",bid_1D,tid_1D,bid,tid_x,tid_y,tid_z); } int main() { const dim3 grid_size(1,65535,1); //大于规定的2^31-1,65535,65535范围虽然能通过编译,但没有输出 const dim3 block_size(1024,1,1); //大于1024虽然能通过编译,但没有输出 hello_from_gpu<<<grid_size,block_size>>>(); cudaDeviceSynchronize(); return 0; }线程束还有一个没有在代码中体现出的概念：线程束，简单理解为同一线程块中相邻的N个线程，一般按架构为32个一组
夫道未始有封，言未始有常，为是而有畛也。作为《深度学习入门》的阅读笔记，本文简略概述感知机，并辅以python实现Python基础再简单的描述离不开基本的语言基础，确保了解《蜻蜓点水-python》中的内容在进行下去逻辑运算和感知机如上所述，python中有逻辑运算，电路中也有逻辑运算，并以“门”的形式呈现，有电流则为1（True），无电流则为0（True）。如下图所示，三种门的符号，真值表就不再赘述，可以用python快速验证感知机的历史1957提出，是神经网络的起源的算法之一。感知机接受多个信号，并输出一个信号。感知机实现先说实现的思路，上述三种逻辑运算都可以用同一种逻辑表达，区别只在于特定参数的不同b称为偏置，w称为权重，二者控制了函数的不同逻辑表示下面直接给出代码实现：初始化函数设定了权重和偏置forward函数在深度学习中代表正向传播，这里借用作为感知机的处理函数逻辑如公式，但将等号给了1class Perceptron: def __init__(self, w=np.array([0,0]),delta=-1.): self.w=w self.delta = delta def forward(self,x): tmp=self.w*x if np.sum(tmp) + self.delta < 0: return 0 else: return 1andclass = Perceptron(w=np.array([0.5,0.5])) print(andclass.forward(np.array([0,0]))) print(andclass.forward(np.array([1,0]))) print(andclass.forward(np.array([0,1]))) print(andclass.forward(np.array([1,1]))) # 0 # 0 # 0 # 1orclass = Perceptron(w=np.array([1,1])) print(orclass.forward(np.array([0,0]))) print(orclass.forward(np.array([1,0]))) print(orclass.forward(np.array([0,1]))) print(orclass.forward(np.array([1,1]))) # 0 # 1 # 1 # 1nandclass = Perceptron(w=np.array([-0.5,-0.5]),delta=0.9) print(nandclass.forward(np.array([0,0]))) print(nandclass.forward(np.array([1,0]))) print(nandclass.forward(np.array([0,1]))) print(nandclass.forward(np.array([1,1]))) # 1 # 1 # 1 # 0目前感知机的局限可以对当前感知机的输入和输出作图，可得：总之，有办法做一条直线将0和1的结果划分开，能这么表示的称为线性感知机但是这样的感知机无法表示XOR操作，也就是异或操作：XOR感知机的实现该问题由1969年被提出，直到1986年被解决，这段时间人工智能的发展是停滞的。解决的方法就是使用多层感知机，即单独的“门”无法解决的问题，就用多个连接一起解决，直接给方法：用代码表示就是:def XOR(x1, x2): s1 = NAND(x1, x2) s2 = OR(x1, x2) y = AND(s1, s2) return y进一步的实现：#def XOR(x1, x2): # s1 = NAND(x1, x2) # s2 = OR(x1, x2) # y = AND(s1, s2) # return y class XOR: def __init__(self): self.nandclass = Perceptron(w=np.array([-0.5,-0.5]),delta=0.9) self.orclass = Perceptron(w=np.array([1,1])) self.andclass = Perceptron(w=np.array([0.5, 0.5])) def forward(self,x): s1 = self.nandclass.forward(x) s2 = self.orclass.forward(x) tmp = self.andclass.forward(np.array([s1,s2])) return tmp print("XOR: ") print(XOR().forward(np.array([0,0]))) print(XOR().forward(np.array([1,0]))) print(XOR().forward(np.array([0,1]))) print(XOR().forward(np.array([1,1]))) #XOR: #0 #1 #1 #0其他问题多层感知机的层数按权重算，2层，因为有两套权重参数按节点（神经元）算，3层，从输入，到中间结果，到输出，有三个阶段多层感知机的进一步感知机通过叠加层能够进行非线性的表示，理论上还可以表示计算机进行的处理多层感知机可以认为是一种神经网络
道在屎溺。作为《深度学习入门》的阅读笔记，本文简略概述python算数计算print(1 - 2) print(4 * 5) print(7 / 5) print(3 ** 2) # 乘方-1 20 1.4 9 数据类型和Type函数type(10)int type(2.718)float type("hello")str type(True)bool 变量x = 100 y = 3.14 print(x * y) type(x * y) # 是注释标识314.0 float 列表a = [1,2,3,4,5] print(a) len(a)[1, 2, 3, 4, 5] 5 print(a[4]) a[4]=99 print(a)5 [1, 2, 3, 4, 99] # 切片 a[0:2] # 获取索引为0到2（不包括2！）的元素[1, 2] a[1:] # 获取从索引为1的元素到最后一个元[2, 3, 4, 99] a[:3] # 获取从第一个元素到索引为3（不包括3！）的元素[1, 2, 3] a[:-1] # 获取从第一个元素到最后一个元素的前一个元素之间的元素[1, 2, 3, 4] a[:-2] # 获取从第一个元素到最后一个元素的前二个元素之间的元素[1, 2, 3] 字典me = {'height':180} # 生成字典 print(me){'height': 180} me['height']=70 print(me['height'])70 布尔型hungry = True sleepy = False type(hungry)bool not hungryFalse hungry and sleepyFalse hungry or sleepyTrue if语句if hungry: print('hungry') elif sleepy: print('sleepy') else: print('not hungry')hungry for 语句for i in [1,2,3]: print(i)1 2 3 函数def hello(object): print("Hello " + object + "!") hello("cat!")Hello cat!! 类class Man: def __init__(self, name): self.name = name print("Initialized!") def hello(self): print("Hello " + self.name + "!") def goodbye(self): print("Good-bye " + self.name + "!") m = Man("David") m.hello() m.goodbye()Initialized! Hello David! Good-bye David! Numpy 速成import numpy as np x = np.array([1.,2.,3.]) print(x) type(x) y = np.array([2.0, 4.0, 6.0])[1. 2. 3.] 以下是对应元素的运算,如果元素个数不同，程序就会报错，所以元素个数保持一致非常重要print(x+y) print(x-y) print(x/y) print(x*y)[3. 6. 9.] [-1. -2. -3.] [0.5 0.5 0.5] [ 2. 8. 18.] NumPy数组不仅可以进行element-wise运算，也可以和单一的数值（标量）组合起来进行运算x/2.0array([0.5, 1. , 1.5]) Numpy的N维数组A = np.array([[1, 2], [3, 4]]) B = np.array([[3, 0],[0, 6]]) print(A)[[1 2] [3 4]] A.shape(2, 2) A.dtypedtype('int32') print(A+B)[[ 4 2] [ 3 10]] print(A*B)[[ 3 0] [ 0 24]] print(A*10)[[10 20] [30 40]] 广播A = np.array([[1, 2], [3, 4]]) B = np.array([10, 20]) A * Barray([[10, 40], [30, 80]]) 访问元素X = np.array([[51, 55], [14, 19], [0, 4]]) X[0][1]55 X[0]array([51, 55]) for row in X: print(row)[51 55] [14 19] [0 4] X = X.flatten() print(X)[51 55 14 19 0 4] 还有一些特殊的访问方法:X[np.array([0, 2, 4])]array([51, 14, 0]) X > 15array([ True, True, False, True, False, False]) X[X>15]array([51, 55, 19]) Matplolib 速成import matplotlib.pyplot as plt # 生成数据 x = np.arange(0, 6, 0.1) # 以0.1为单位，生成0到6的数据 y1 = np.sin(x) y2 = np.cos(x) # 绘制图形 plt.plot(x, y1) plt.show()​​plt.plot(x, y1, label="sin") plt.plot(x, y2, linestyle = "--", label="cos") # 用虚线绘制 plt.xlabel("x") # x轴标签 plt.ylabel("y") # y轴标签 plt.title('sin & cos') # 标题 plt.legend() #图例 plt.show()​ ​from matplotlib.image import imread img = imread('../dataset/lena.png') # 读入图像（设定合适的路径！） plt.imshow(img) plt.show()​ ​
原因在调整spark配置问题时把集群内部通信接口暴露给了公网，提供了攻击的机会攻击者通过端口扫描的方式，用暴露的接口【通过某种方式】进入集群，集群内部可无密码ssh访问，直接访问ssh目录生成密钥，并修改authorized_keys权限添加redis和程序包（注入+进程隐藏）https://github.com/arget13/memdlopen/blob/main/README.mdhttps://github.com/gianlucaborello/libprocesshider/blob/master/README.md这导致了top无法发现其真正进程发现cpu占用异常，单核100%，和top无法对应netstat -alntp 发现异常tcp链接：5.180.174.211:80，检查为境外乌克兰ip，lastb发现有大量异常失败登录，且ip大多来自境外，疑似撞库处理关掉集群其他机器立即修改安全组，关闭除22以外的公网端口试图删除authorized_keys known_host中的攻击者数据，发现并无权限，但查看权限并无异常，使用lsattr查看文件属性发现被添加ia标记，无法修改无法删除，使用chattr 修改无效且返回结果看起来比较诡异，查看/usr/bin/chattr 发现已被修改，使用ls -l |grep Jan 过滤出最近修改的所有文件，发现其他新增命令工具但此时需要chattr 来修改标识才能进行操作，从https://github.com/posborne/linux-programming-interface-exercises/blob/master/15-file-attributes/chattr.c下载文件gcc编译，重命名为chattr，利用其来修改标识 chattr -ia，删除原本的chattr 替换之，然后就可以对其他文件进行操作https://c.biancheng.net/view/874.html查看/etc/ld.so.preload文件，此文件原先应为空或不存在，设法清空或删除该文件（和文件里的文件路径），top的问题就可以恢复参考：https://www.jianshu.com/p/d4b3d87fdc51也可通过busybox top工具直接发现真正的进程，定位到程序所在的目录（我的是/var/tmp/.11/，里面都是编译后的二进制文件和上述github的bash脚本），先不急删使用systemctl status 对应pid，发现有对应的service（于我情况如此，也可能会使用定时程序），清理掉就可以，然后附带检查service而文件的修改日期，发现有其他攻击者新增但未启动的service一并处理另外本环境没有redis服务，故连带清理掉，这样cpu和内存就应该能降下来最后删除病毒文件即可然后关闭该机器，对其他机器依次操作（因为集群是免密互通的，保险起见）后续改进关闭所有不必要的公网端口内部web服务使用nginx反代并管理登录Nginx 配置 安全认证 反向代理 HDFS web 页面_nginx如何代理一个web页面-CSDN博客