"iteration over a 0-d array" error in multiprocessing when applying zip function

Question

Dear programmers: We tried to solve a non-linear coupled Poisson drift-diffusion PDE in Fipy, which implemented well. But, considering the amount of processing(which takes less than a month) we tried to solve the equations in parallel via multiprocessing, therefore we made a file for solving equations(EQN) and import it into the main file in jupyterLab. Having four variables we used the zip function, but despite converting the input of starmap to list, it gives ” iteration over a 0-d array” error which I think is due to CellVariables but not sure how to solve it. Any help regarding this problem or tip regarding otherways of multiprocessing is highly appreciated.

The main code is summarized as follows:

import EQN
from itertools import zip_longest
from fipy import *
from fipy import Grid1D, CellVariable, TransientTerm, DiffusionTerm, Viewer
from fipy.solvers import *
from multiprocessing import Pool
import numpy as np
import math
import matplotlib.pyplot as plt
from matplotlib import cm
FIPY_SOLVERS=scipy
#____________________________________________________
#Defining main constants
 ……….
#____________________________________________________
#Defining Variables
m = Grid1D(nx=nx, dx=dx, overlap=2, communicator=parallelComm)
phi = CellVariable(mesh=m, hasOld=True, value=0.)
ne = CellVariable(mesh=m, hasOld=True, value=1.0)
phi_face = phi.faceValue
ne_face = ne.faceValue
x = m.cellCenters[0]
t0 = Variable()
phi.setValue(0)
ne.setValue(1.0)
#__________________________
#Boundary Condition:
valueleft_phi= 0
valueright_phi= 0
valueleft_ne=1.0
valueright_ne= 1.0
phi.constrain(valueleft_phi, m.facesLeft)
phi.constrain(valueleft_phi, m.facesRight)
ne.constrain(valueright_ne, m.facesLeft)
ne.constrain(valueright_ne, m.facesRight)
#____________________________________________________
steps = 3e7
T=4
dt=T/steps
F=dt/(dx**2)
print('F=',F)
#____________________________________________________
#Main Line
vi = Viewer(ne)
with open('out2.txt', 'w') as output:
    while t0()<T:
        print(t0)
        phi.updateOld()
        ne.updateOld()
        res=1e40
        while res > 1e5:
            if __name__=='__main__':
                task=list(zip_longest(phi,ne,x,t0))
                num_processors = 2
                p = Pool(processes=num_processors)
                res=p.starmap(EQN.PHI_NE,task)
                mySolver = LinearGMRESSolver
            t0.setValue(t0() + dt)

Below is the summarized python function code(EQN) which we imported into the main JupyterLab file:

from fipy import *
from fipy import Grid1D, CellVariable, TransientTerm, DiffusionTerm
from fipy.solvers import *
import numpy as np
import math
#----------------------------------------------------------------------
#Defining Constants
  ………
#--------------------------------------------------------------
def S(x,t):
    f=-(n0/n_i)*(e**(-(np.log(2))*(((x-d0)/X0)**2)))*(1-e**(-t))*np.heaviside((1-t),0)
    return f
def PHI_NE(phi,ne,x,t):
    eqn0 =  DiffusionTerm(1.,var=phi)==ImplicitSourceTerm(1,var=ne)-1.0
    eqn1 = TransientTerm(1.,var=ne) ==  - VanLeerConvectionTerm(A*phi.faceGrad\
                                                       ,var=ne)+DiffusionTerm(B,var=ne)+S(x,t)

    eqn = eqn0 & eqn1
    return eqn.sweep