def step(self, x, states):
		M = states[0]  # (nb_samples, nb_slots, memory_size)
		h = states[1]  # (nb_samples, memory_size)
		w = states[2]  # (nb_samples, nb_slots)
		#------Memory read--------#
		k = self.W_k(h)  # (nb_samples, memory_size)
		w_hat = T.batched_tensordot(M, k, axes=[(2), (1)])  # (nb_samples, nb_slots)
		beta = K.sigmoid(self.W_b(h))  # (nb_samples, 1)
		beta = K.repeat(beta, self.nb_slots)  # (nb_samples, nb_slots, 1)
		beta = K.squeeze(beta, 2)  # (nb_samples, nb_slots)
		w_hat = softmax(w_hat * beta)  # (nb_samples, nb_slots)
		g = sigmoid(self.W_hg(h))  # (nb_samples, 1)
		g = K.repeat(g, self.nb_slots)  # (nb_samples, nb_slots, 1)
		g = K.squeeze(g, 2)  # (nb_samples, nb_slots)
		w = (1 - g) * w + g * w_hat  # (nb_samples, nb_slots)
		c = T.batched_tensordot(w, M, axes=[(1), (1)])
		h = tanh(self.W_ih(x) + self.W_c(c))
		y = self.W_ho(h)
		#---------Memory write---------#
		v = self.W_v(h)  # (nb_samples, memory_size)
		v = K.repeat(v, 1)
		e = sigmoid(self.W_he(h))  # (nb_samples, nb_slots)
		f = 1 - w * e  # (nb_samples, nb_slots)
		f = K.repeat(f, self.memory_size)  # (nb_samples, memory_size, nb_slots)
		f = K.permute_dimensions(f, (0, 2, 1))  # (nb_samples, nb_slots, memory_size)
		u = w  # (nb_samples, nb_slots)
		u = K.repeat(u, 1)
		uv = T.batched_tensordot(u, v, axes=[(1), (1)])
		M = M * f + uv
		return y, [M, h, w]

    def __init__(self, intpic_parameters=None,
            case_costs=None, pics=None, case_labels=None,
            batch_size=None, pic_size=None, label_count=None, **kwargs):
        super(IntpicGradientDescent, self).__init__(**kwargs)
        center_val = 0.5
        self.input_pics = pics
        self.case_costs = case_costs
        self.batch_size = batch_size
        self.label_count = label_count
        self.intpic_parameters = intpic_parameters
        self.jacobians = self._compute_jacobians()
        self.gradpics = OrderedDict(
          [(param, _create_intpic_histogram_for(param, pic_size, label_count))
                for param in self.intpic_parameters])
        self.intpics = OrderedDict(
          [(param, _create_intpic_histogram_for(param, pic_size, label_count))
                for param in self.intpic_parameters])
        # attributes pics: (cases, picy, picx) to (cases, labels, picy, picx)
        # attributed_pics = tensor.batched_tensordot(
        #     tensor.extra_ops.to_one_hot(case_labels.flatten(), label_count),
        #     pics[:, 0, :, :], axes=0)
        zeroed_pics = pics - 0.5
        attributed_pics = tensor.batched_tensordot(
            tensor.extra_ops.to_one_hot(
                case_labels.flatten(), label_count),
            zeroed_pics[:, 0, :, :],
            axes=0)
        self.gradpic_updates = OrderedDict(
            [_create_gradpic_updates(
                self.gradpics[param],
                self.jacobians[param],
                attributed_pics) for param in self.intpic_parameters])
        self.add_updates(self.gradpic_updates)

        intensity_pics = (zeroed_pics *
                gradient.grad(case_costs.mean(), pics))
        attributed_i_pics = tensor.batched_tensordot(
            tensor.extra_ops.to_one_hot(
                case_labels.flatten(), label_count),
            intensity_pics[:, 0, :, :],
            axes=0)

        self.intpic_updates = OrderedDict(
            [_create_intensity_updates(
                self.intpics[param],
                self.jacobians[param],
                attributed_i_pics) for param in self.intpic_parameters])
        self.add_updates(self.intpic_updates)

    def negLeftFactorization(self, batchSize, negEmbA, argsEmbB, wC, wC1, wC2):
        # l = batchSize
        # k = self.k  # embed size
        # r = self.r  # relation number
        # argEmbedsA = self.A[argsA.flatten()]  # [l,k]
        # argEmbedsB = self.A[argsB.flatten()]  # [l,k]
        # first = T.tensordot(relationProbs, self.C, axes=[[1], [2]])  # [l,r] * [k,k,r] = [l, k, k]
        Afirst = T.batched_tensordot(wC, negEmbA.dimshuffle(1, 2, 0), axes=[[1], [1]])  # [l, k, k] * [n, l, k] = [l, k, n]
        Asecond = T.batched_tensordot(Afirst, argsEmbB, axes=[[1], [1]])  # [l, k, n] * [l, k] = [l, n]

        # spFirst = T.dot(relationProbs, self.C1.dimshuffle(1, 0))  # [l,r] * [r,k] = [l, k]
        spAfirst = T.batched_tensordot(wC1, negEmbA.dimshuffle(1, 2, 0), axes=[[1], [1]])  # [l,k] [l,k,n] = [l,n]

        spSecond = T.batched_dot(wC2, argsEmbB)

        return Asecond + spAfirst + spSecond.reshape((batchSize, 1))

    def get_output_for(self, inputs, **flags):
        u, mu, L = inputs

        # batch-wise matrix multiplication P = L * L.T
        P = T.batched_tensordot(L, L.swapaxes(2, 1), axes=[2, 1])

        # (u-mu) * P for each batch
        diff_times_P = T.batched_tensordot((u - mu), P, axes=[1, 2])

        # A = - 0.5 * (u-mu) * P * (u-mu).T
        A = -0.5 * T.batched_dot(diff_times_P, (u - mu))[:,None]

        #shape = (None,1)
        assert A.ndim ==2

        return A

def batch_dot(x, y, axes=None):
    '''batchwise dot product
    batch_dot results in a tensor with less dimensions than the input.
    If the number of dimensions is reduced to 1, we use `expand_dims` to
    make sure that ndim is at least 2.

    # Example
        Assume x = [[1, 2], [3, 4]]   and y = [[5, 6], [7, 8]]
        batch_dot(x, y, axes=1) = [[17, 53]] which is the main diagonal
        of x.dot(y.T), although we never have to calculate the off-diagonal
        elements.


    # Arguments
        x, y: tensors with ndim >= 2
        axes: list (or single) int with target dimensions

    # Returns
        Tensor with ndim >= 2
    '''
    if type(axes) == int:
        axes = (axes, axes)
    if axes is None:
        # behaves like tf.batch_matmul as default
        axes = [x.ndim - 1, y.ndim - 2]
    out = T.batched_tensordot(x, y, axes=axes)
    if ndim(out) == 1:
        out = expand_dims(out, 1)
    return out

 def negLeftMostFactorization(self, batchSize, negEmbed, wC1):
     # l = batchSize
     # k = self.k  # embed size
     # r = self.r  # relation number
     # first = T.dot(relationProbs, self.C1.dimshuffle(1, 0))  # [l,r] * [r,k] = [l, k]
     Afirst = T.batched_tensordot(wC1, negEmbed.dimshuffle(1, 2, 0), axes=[[1], [1]])  # [l,k] [l,k,n] = [l,n]
     return Afirst

    def factorization(self, batchSize, argsEmbA, argsEmbB, wC):

        # first = T.tensordot(relationProbs, self.C, axes=[[1], [2]])  # [l,r] * [k,k,r] = [l, k, k]
        Afirst = T.batched_tensordot(wC, argsEmbA, axes=[[1], [1]])  # [l, k, k] * [l, k] = [l, k]
        Asecond = T.batched_dot(Afirst, argsEmbB)  # [l, k] * [l, k] = [l]
        # entropy = T.sum(T.log(relationProbs) * relationProbs, axis=1)  # [l,r] * [l,r] = [l]
        return Asecond

 def negRightMostFactorization(self, batchSize, negEmbed, wC2):
     # l = batchSize
     # k = self.k  # embed size
     # r = self.r  # relation number
     # second = T.dot(relationProbs, self.C2.dimshuffle(1, 0))  # [l,r] * [r,k] = [l, k]
     Asecond = T.batched_tensordot(wC2, negEmbed.dimshuffle(1, 2, 0), axes=[[1], [1]])  # [l,k] [l,k,n] = [l,n]
     return Asecond

 def attention(self, m, q, mask):
     # mask original shape is (batch*memory_length, input_length, 1)
     # shape (batch, memory)
     mask = K.reshape(mask[:, 0], (-1, self.memory_length))
     # shape: (batch, memory_length, 1)
     p = T.batched_tensordot(m, q, (2, 2))
     # shape: (batch, memory_length)
     p = K.softmax(p[:, :, 0])  # * K.cast(mask, 'float32')
     # shape: (batch, 1, memory_length)
     return K.expand_dims(p, dim=1)

def batch_matmul(x, y, adj_x = False, adj_y = False):
  # used for tensor, i.e., if x and y are matrix, then use matrix
  # multiplication, else batch_matmul
  # TODO: interfaces with tensorflow backed
  axes = (x.ndim - 1, y.ndim - 2)
  if adj_x == True:
    axes[0] = x.ndim - 2
  if adj_y == True:
    axes[1] = y.ndim - 1
  return T.batched_tensordot(x, y, axes = axes)

    def get_output_for(self, inputs, **kwargs):
        """

        :param inputs: inputs: list of theano.TensorType
            `inputs[0]` should always be the symbolic input variable.  When
            this layer has a mask input (i.e. was instantiated with
            `mask_input != None`, indicating that the lengths of sequences in
            each batch vary), `inputs` should have length 2, where `inputs[1]`
            is the `mask`.  The `mask` should be supplied as a Theano variable
            denoting whether each time step in each sequence in the batch is
            part of the sequence or not.  `mask` should be a matrix of shape
            ``(n_batch, n_time_steps)`` where ``mask[i, j] = 1`` when ``j <=
            (length of sequence i)`` and ``mask[i, j] = 0`` when ``j > (length
            of sequence i)``.
        :return: theano.TensorType
            Symbolic output variable.
        """
        input = inputs[0]
        mask = None
        if self.mask_incoming_index > 0:
            mask = inputs[self.mask_incoming_index]

        # compute the bi-affine part
        # first via tensor dot ([batch, length, dim] * [dim, dim, num_label])
        # output shape = [batch, length, dim, num_label]
        out = T.tensordot(input, self.U, axes=[[2], [0]])
        # second via tensor dot ([batch, length, dim, num_label] * [batch, dim, length)
        # output shape = [batch, length, length, num_label]
        out = T.batched_tensordot(out, input.dimshuffle(0, 2, 1), axes=([2], [1]))
        out = out.dimshuffle(0, 1, 3, 2)

        # compute head bias part by tensor dot ([batch, length, dim] * [dim, num_label])
        # the shape of s_h should be [batch, length, num_label]
        if self.W_h is not None:
            s_h = T.tensordot(input, self.W_h, axes=[[2], [0]])
            out = out + s_h.dimshuffle(0, 1, 'x', 2)

        # compute child part by tensor dot ([batch, length, dim] * [dim, num_label]
        # the shape of s_c should be [batch, length, num_label]
        if self.W_c is not None:
            s_c = T.tensordot(input, self.W_c, axes=[[2], [0]])
            out = out + s_c.dimshuffle(0, 'x', 1, 2)

        # add bias part.
        if self.b is not None:
            out = out + self.b.dimshuffle('x', 'x', 'x', 0)

        if mask is not None:
            mask_shuffled = mask.dimshuffle(0, 1, 'x', 'x')
            out = out * mask_shuffled
            mask_shuffled = mask.dimshuffle(0, 'x', 1, 'x')
            out = out * mask_shuffled
        return out

    def factorization(self, batchSize, argsEmbA, argsEmbB, wC, wC1, wC2):
        # l = batchSize
        # k = self.k  # embed size
        # r = self.r  # relation number
        # argEmbedsA = self.A[argsA.flatten()]  # [l,k]
        # argEmbedsB = self.A[argsB.flatten()]  # [l,k]

        # first = T.tensordot(relationProbs, self.C, axes=[[1], [2]])  # [l,r] * [k,k,r] = [l, k, k]
        Afirst = T.batched_tensordot(wC, argsEmbA, axes=[[1], [1]])  # + self.Cb  # [l, k, k] * [l, k] = [l, k]
        Asecond = T.batched_dot(Afirst, argsEmbB)  # [l, k] * [l, k] = [l]
        # entropy = T.sum(T.log(relationProbs) * relationProbs, axis=1)  # [l,r] * [l,r] = [l]
        spFirst = T.batched_dot(wC1, argsEmbA)
        spSecond = T.batched_dot(wC2, argsEmbB)
        return Asecond + spFirst + spSecond

 def output(self, input_vectors, input_scalars):
     """
     Calculate the n_output transformed vectors for this layer
     @param input_scalars: n_input x n_output scalar vector
     @param input_vectors: n_input vectors (actual shape should be (n_batch, n_input, n_dimension)
     """
     mat = input_scalars.reshape((n_batch, self.n_input, self.n_output))
     z = T.batched_tensordot(input_vectors, mat, [[1], [1]]).swapaxes(1, 2) + T.addbroadcast(self.b, 0, 2)
     if self.activation == 'linear':
         return z
     elif self.activation == 'rectified':
         return T.maximum(z, 0)
     elif self.activation == 'tanh':
         return T.tanh(z)
     else:
         raise "Unknown activation, %s" % self.activation

def batch_dot(x, y, axes=None):
    '''Batchwise dot product.

    batch_dot results in a tensor with less dimensions than the input.
    If the number of dimensions is reduced to 1, we use `expand_dims` to
    make sure that ndim is at least 2.

    # Arguments
        x, y: tensors with ndim >= 2
        axes: list (or single) int with target dimensions

    # Returns
        A tensor with shape equal to the concatenation of x's shape
        (less the dimension that was summed over) and y's shape
        (less the batch dimension and the dimension that was summed over).
        If the final rank is 1, we reshape it to (batch_size, 1).

    # Examples
        Assume x = [[1, 2], [3, 4]]   and y = [[5, 6], [7, 8]]
        batch_dot(x, y, axes=1) = [[17, 53]] which is the main diagonal
        of x.dot(y.T), although we never have to calculate the off-diagonal
        elements.

        Shape inference:
        Let x's shape be (100, 20) and y's shape be (100, 30, 20).
        If dot_axes is (1, 2), to find the output shape of resultant tensor,
            loop through each dimension in x's shape and y's shape:
        x.shape[0] : 100 : append to output shape
        x.shape[1] : 20 : do not append to output shape,
            dimension 1 of x has been summed over. (dot_axes[0] = 1)
        y.shape[0] : 100 : do not append to output shape,
            always ignore first dimension of y
        y.shape[1] : 30 : append to output shape
        y.shape[2] : 20 : do not append to output shape,
            dimension 2 of y has been summed over. (dot_axes[1] = 2)

        output_shape = (100, 30)
    '''
    if type(axes) == int:
        axes = (axes, axes)
    if axes is None:
        # behaves like tf.batch_matmul as default
        axes = [x.ndim - 1, y.ndim - 2]
    out = T.batched_tensordot(x, y, axes=axes)
    if ndim(out) == 1:
        out = expand_dims(out, 1)
    return out

 def get_output(self, train=False):
     if self.mode == 'sum' or self.mode == 'ave':
         s = self.layers[0].get_output(train)
         for i in range(1, len(self.layers)):
             s += self.layers[i].get_output(train)
         if self.mode == 'ave':
             s /= len(self.layers)
         return s
     elif self.mode == 'concat':
         inputs = [self.layers[i].get_output(train) for i in range(len(self.layers))]
         return T.concatenate(inputs, axis=self.concat_axis)
     elif self.mode == 'join':
         inputs = OrderedDict()
         for i in range(len(self.layers)):
             X = self.layers[i].get_output(train)
             if X.name is None:
                 raise ValueError("merge_mode='join' only works with named inputs")
             else:
                 inputs[X.name] = X
         return inputs
     elif self.mode == 'mul':
         s = self.layers[0].get_output(train)
         for i in range(1, len(self.layers)):
             s *= self.layers[i].get_output(train)
         return s
     elif self.mode == 'dot':
         l1 = self.layers[0].get_output(train)
         l2 = self.layers[1].get_output(train)
         output = T.batched_tensordot(l1, l2, self.dot_axes)
         output_shape = list(self.output_shape)
         output_shape[0] = l1.shape[0]
         output = output.reshape(tuple(output_shape))
         return output
     elif self.mode == 'inner':
         l1 = self.layers[0].get_output(train)
         l2 = self.layers[1].get_output(train)
         output =T.sum(l1 * l2, axis=-1)
         return output
     elif self.mode == 'cos':
         l1 = self.layers[0].get_output(train)
         l2 = self.layers[1].get_output(train)
         output, _ = theano.scan(lambda v1, v2: T.dot(v1, v2) / T.sqrt(T.dot(v1, v1) * T.dot(v2, v2)),
                                 sequences=[l1, l2],
                                 outputs_info=None)
         return output
     else:
         raise Exception('Unknown merge mode')

def weighted_average(inp, weights, axis=None):
    # n_b x n_s x 4 x n_w_a: inp
    if axis == 2:  # for question
        weights = weights.flatten(ndim=2)
        weights /= T.sum(weights, axis=1, keepdims=True) + 0.000001
        return T.batched_tensordot(inp, weights, [[inp.ndim - 1], [1]])
    elif axis == 3:  # for answer inp: (None, 51, 4, 20), output: (None, 4, 20, 1)
        weights = weights.flatten(ndim=weights.ndim - 1)
        weights /= T.sum(weights, axis=weights.ndim - 1, keepdims=True) + 0.000001
        weights = weights.dimshuffle(0, 'x', 1, 2)
        return T.sum(inp * weights, axis=3)
    elif axis == 4:  # for inner sliding window
        weights = weights.flatten(ndim=weights.ndim - 1)
        weights /= T.sum(weights, axis=weights.ndim - 1, keepdims=True) + 0.000001
        weights = weights.dimshuffle(0, 'x', 'x', 1, 2)
        return T.sum(inp * weights, axis=4)
    else:
        raise RuntimeError

 def get_output(self, train=False):
     if self.mode == "sum" or self.mode == "ave":
         s = self.layers[0].get_output(train)
         for i in range(1, len(self.layers)):
             s += self.layers[i].get_output(train)
         if self.mode == "ave":
             s /= len(self.layers)
         return s
     elif self.mode == "concat":
         inputs = [self.layers[i].get_output(train) for i in range(len(self.layers))]
         return T.concatenate(inputs, axis=self.concat_axis)
     elif self.mode == "join":
         inputs = OrderedDict()
         for i in range(len(self.layers)):
             X = self.layers[i].get_output(train)
             if X.name is None:
                 raise ValueError("merge_mode='join' only works with named inputs")
             else:
                 inputs[X.name] = X
         return inputs
     elif self.mode == "mul":
         s = self.layers[0].get_output(train)
         for i in range(1, len(self.layers)):
             s *= self.layers[i].get_output(train)
         return s
     elif self.mode == "dot":
         l1 = self.layers[0].get_output(train)
         l2 = self.layers[1].get_output(train)
         output = T.batched_tensordot(l1, l2, self.dot_axes)
         return output
     elif self.mode == "cos":
         l1 = self.layers[0].get_output(train)
         l2 = self.layers[1].get_output(train)
         output, _ = theano.scan(
             lambda v1, v2: T.dot(v1, v2) / T.sqrt(T.dot(v1, v1) * T.dot(v2, v2)),
             sequences=[l1, l2],
             outputs_info=None,
         )
         return output
     else:
         raise Exception("Unknown merge mode")

        def step(i, in_mask, ACT, ACT_, in_se, WT):
            sub_tree_idx_ = T.nonzero(WT[:, i, :] > -1)
            a_ = T.dot(in_se[:, i], self.WSM)  # + self.b
            if self.b is not None:
                a_ += self.b.dimshuffle('x', 0)
            a_ = a_ + T.sum(ACT_[:, i], axis=1)
            a_ = T.tanh(a_)
#            if self.dropout:
#                a_ = a_ / self.retain_prob * self._srng.binomial(a_.shape, p=self.retain_prob,
#                                                                 dtype=theano.config.floatX)
            a_ = T.switch(in_mask, a_, ACT[:, i-1])
            a__ = T.batched_tensordot(a_[sub_tree_idx_[0], :],
                                      self.WC[WT[sub_tree_idx_[0], 
                                                 i, sub_tree_idx_[1]]], axes=1)
#            if self.dropout:
#                a__ = a__ / self.retain_prob * self._srng.binomial(a__.shape, p=self.retain_prob,
#                                                                   dtype=theano.config.floatX)
            newACT_ = T.set_subtensor(ACT_[sub_tree_idx_[0], sub_tree_idx_[1], i],
                                      a__)
            newACT = T.set_subtensor(ACT[:, i], a_)
            return newACT, newACT_

 def __init__(self, actpic_variables=None, pics=None, case_labels=None,
         label_count=None, data_stream=None, rectify=False, **kwargs):
     center_val = 0.5
     self.input_pics = pics
     # self.batch_size = batch_size
     # self.label_count = label_count
     self.actpic_variables = actpic_variables
     # attributes pics: (cases, picy, picx) to (cases, labels, picy, picx)
     # attributed_pics = tensor.batched_tensordot(
     #     tensor.extra_ops.to_one_hot(case_labels.flatten(), label_count),
     #     pics[:, 0, :, :], axes=0)
     zeroed_pics = pics - 0.5
     attributed_pics = tensor.batched_tensordot(
         tensor.extra_ops.to_one_hot(
             case_labels.flatten(), label_count),
         zeroed_pics[:, 0, :, :],
         axes=0)
     self.actpics = [self._create_actpic_image_for(
           name + '_actpic', var, attributed_pics, rectify)
             for name, var in self.actpic_variables.items()]
     self.evaluator = DatasetEvaluator(self.actpics)
     self.data_stream = data_stream
     self.results = None
     super(ActpicExtension, self).__init__(**kwargs)

    def get_output(self, train=False):
        if self.mode == 'sum' or self.mode == 'ave':
            s = self.layers[0].get_output(train)
            for i in range(1, len(self.layers)):
                s += self.layers[i].get_output(train)
            if self.mode == 'ave':
                s /= len(self.layers)
            return s
        elif self.mode == 'index':
            return self.layers[0].get_output(train)[self.layers[1].get_output(train)[:, 0]]
        elif self.mode == 'imax' or self.mode == 'iavg' or self.mode == 'imaxavg' or \
                self.mode =='imaxminavg':
            if self.mode == 'imax':
                fn = lambda start, end, x: T.max(x[start[0]:end[0]], axis=0)
            elif self.mode == 'iavg':
                fn = lambda start, end, x: T.mean(x[start[0]:end[0]], axis=0)
            elif self.mode == 'imaxavg':
                fn = lambda start, end, x: T.concatenate([T.max(x[start[0]:end[0]], axis=0),
                                                          T.mean(x[start[0]:end[0]], axis=0)])
            else:
                fn = lambda start, end, x: T.concatenate([T.max(x[start[0]:end[0]], axis=0),
                                                          T.min(x[start[0]:end[0]], axis=0),
                                                          T.mean(x[start[0]:end[0]], axis=0)])

            data = self.layers[0].get_output(train)
            starts = self.layers[1].get_output(train)
            ends = self.layers[2].get_output(train)
            outputs, _ = theano.scan(fn=fn,
                                     outputs_info=None,
                                     sequences=[starts, ends],
                                     non_sequences=data)
            return outputs
        elif self.mode == 'mm':
            scores = self.layers[0].get_output(train)
            starts = self.layers[1].get_output(train)
            ends = self.layers[2].get_output(train)
            costs = self.layers[3].get_output(train)
            fn = lambda start, end, scs, csts: \
                T.max(csts[start[0]:end[0]] *
                      (3 + scs[start[0]:end[0]] -
                       T.max(scs[start[0]:end[0]][T.eq(csts[start[0]:end[0]], 0).nonzero()[0]])))
            outputs, _ = theano.scan(fn=fn,
                                     outputs_info=None,
                                     sequences=[starts, ends],
                                     non_sequences=[scores, costs])
            return outputs.reshape((outputs.size, 1))
        elif self.mode == 'risk':
            scores = self.layers[0].get_output(train)
            starts = self.layers[1].get_output(train)
            ends = self.layers[2].get_output(train)
            costs = self.layers[3].get_output(train)
            fn = lambda start, end, scs, csts: \
                T.sum(costs[start[0]:end[0]] * T.nnet.softmax(scs[start[0]:end[0]].T).T)
            outputs, _ = theano.scan(fn=fn,
                                     sequences=[starts, ends],
                                     non_sequences=[scores, costs])
            return outputs.reshape((outputs.size, 1))
        elif self.mode == 'concat':
            inputs = [self.layers[i].get_output(train) for i in range(len(self.layers))]
            return T.concatenate(inputs, axis=self.concat_axis)
        elif self.mode == 'join':
            inputs = OrderedDict()
            for i in range(len(self.layers)):
                X = self.layers[i].get_output(train)
                if X.name is None:
                    raise ValueError("merge_mode='join' only works with named inputs")
                else:
                    inputs[X.name] = X
            return inputs
        elif self.mode == 'mul':
            s = self.layers[0].get_output(train)
            for i in range(1, len(self.layers)):
                s *= self.layers[i].get_output(train)
            return s
        elif self.mode == 'dot':
            l1 = self.layers[0].get_output(train)
            l2 = self.layers[1].get_output(train)
            output = T.batched_tensordot(l1, l2, self.dot_axes)
            output_shape = list(self.output_shape)
            output_shape[0] = l1.shape[0]
            output = output.reshape(tuple(output_shape))
            return output
        elif self.mode == 'cos':
            l1 = self.layers[0].get_output(train)
            l2 = self.layers[1].get_output(train)
            output, _ = theano.scan(lambda v1, v2: T.dot(v1, v2) / T.sqrt(T.dot(v1, v1) * T.dot(v2, v2)),
                                    sequences=[l1, l2],
                                    outputs_info=None)
            return output
        else:
            raise Exception('Unknown merge mode')