So perhaps something like this?
def constraint(y,param_mat):
    return torch.mm(param_mat,y)

param_mat = torch.Parameter(torch.tensor(np.random.normal(0,1,(y_size,y_size)),requires_grad = True)
param_list = [list(model_list[i].parameters()) + list(params_mat) for i in range(num_models)]
optimizers_list = [optim.Adam(param_list[i]) for i in range(num_models)]
for i in range(num_nets):
    ipt = batch[0][i]
    target = batch[1][i]
    y = model_list[i](ipt)
    loss = loss1(y,target)+loss2(constraint(t,param_mat))
    optimizers_list[i].zero_grad()
    loss.backwards()
    optimizers_list[i].step()