code: def update_gradients_full(self, dL_dK, X, X2=None):
          
          if X2 is None:
              X2 = X
  
          base = np.pi * (X[:, None, :] - X2[None, :, :]) / self.period
  
          sin_base = np.sin( base )
          exp_dist = np.exp( -0.5* np.sum( np.square(  sin_base / self.lengthscale ), axis = -1 ) )
  
          dwl = self.variance * (1.0/np.square(self.lengthscale)) * sin_base*np.cos(base) * (base / self.period)
  
          dl = self.variance * np.square( sin_base) / np.power( self.lengthscale, 3)
  
          self.variance.gradient = np.sum(exp_dist * dL_dK)
          #target[0] += np.sum( exp_dist * dL_dK)
  
          if self.ARD1: # different periods
              self.period.gradient = (dwl * exp_dist[:,:,None] * dL_dK[:, :, None]).sum(0).sum(0)
          else:  # same period
              self.period.gradient = np.sum(dwl.sum(-1) * exp_dist * dL_dK)
  
          if self.ARD2: # different lengthscales
              self.lengthscale.gradient = (dl * exp_dist[:,:,None] * dL_dK[:, :, None]).sum(0).sum(0)
          else: # same lengthscales
              self.lengthscale.gradient = np.sum(dl.sum(-1) * exp_dist * dL_dK)
docs: derivative of the covariance matrix with respect to the parameters.