As far as I know, I am not changing the dimension of toep. The dimensions of this list in my case is shown below which I observed to be uniform across all input data.
toep.len(): 30
toep[0].shape: torch.Size([12, 2, 320, 320])
I am not aware if I am consciously storing references to this variable so that its deletion is prevented. However, I am fairly confident that after taking this variable back to CPU as shown above, it really solves the problem (or workarounds)
This is my workflow
Our operator which we apply in the CG iterations depends on the data. This operator is computed just before we start our training process. But because we are limited with our GPU memory, we store components of this operator in CPU memory. We only call it in GPU for the respective batch of data when we run our network.
# Creating our operator 
for i in range(total_data): 
   csm = load_coil_map(csmsfolder, csmfname)
   csm = torch.tensor(csm).cuda()       ​
   ​enc_obj = Dyn2DRadEncObj(csm , .. other parameters ..)
   ​for j in range(len(enc_obj.AdagA_toep_kernel_list)):
       ​enc_obj.AdagA_toep_kernel_list[j] = enc_obj.AdagA_toep_kernel_list[j].cpu()
   ​op[i] = [csm_tensor.cpu(), enc_obj.AdagA_toep_kernel_list]

Next, we call our net during our training process. Notice we set the opcomp attribute of the net with the respective op index.
with torch.no_grad():
    ​for index in range(Ids): 
       bx = x[index,:]    # the input
       by = y[index,:]    # the labels
       bop = op[index,:]  # data dependent operator
       ​net.opcomp = bop
       ​l = loss(net(bx.detach.()),by)
       ​.
       ​.
   ​.
   ​.
​.

Inside the network we assign it to csm and toep just before doing the conjugate gradient. After that, we do a backpropagation to learn the netwok weights.
   ​for mb_indx in range(mini_batch_size):
       csm = self.opcomp[0].cuda()
   ​    toep = self.opcomp[1]
       .
       xm = self.ConjGrad(xm, ..other parameters..)

.

I am sorry if it sounds like I am complicating things but our problem is that the CG function at the moment only works with mini-batch = 1. For mini batchsize more that one, we call a loop for the conjugate gradient, and self.opcomp is the list of operator components corresponding to the required batch.