def METHOD_NAME(self, model):
    print('=' * 100)
    if self.complex():
        modelRe = model[0:int(len(model)/2)]
        modelIm = model[int(len(model)/2):len(model)]
        modelC = pg.math.toComplex(modelRe, modelIm)
        print("Real", min(modelRe), max(modelRe))
        print("Imag", min(modelIm), max(modelIm))
        u = self.prepareJacobian_(modelC)
        if self._J.rows() == 0:
            #re(data)/re(mod) = im(data)/im(mod)
            # we need a local copy until we have a gimli internal reference counter FIXTHIS
            M1 = pg.Matrix()
            M2 = pg.Matrix()
            self.matrixHeap.append(M1)
            self.matrixHeap.append(M2)
            JRe = self._J.addMatrix(M1)
            JIm = self._J.addMatrix(M2)
            self._J.addMatrixEntry(JRe, 0, 0)
            self._J.addMatrixEntry(JIm, 0, len(modelRe), -1.0)
            self._J.addMatrixEntry(JIm, self.data().size(), 0, 1.0)
            self._J.addMatrixEntry(JRe, self.data().size(), len(modelRe))
        else:
            self._J.clean()
        k = pg.Vector(self.data()('k'))
        self.data().set('k', k*0.0 + 1.0)
        dMapResponse = pb.DataMap()
        dMapResponse.collect(self.electrodes(), self.solution())
        respRe = dMapResponse.data(self.data(), False, False)
        respIm = dMapResponse.data(self.data(), False, True)
        #CVector resp(toComplex(respRe, respIm));
        #RVector am(abs(resp) * dataContainer_->get("k"));
        #RVector ph(-phase(resp));
        print("respRe", pg.math.median(respRe), min(respRe), max(respRe))
        print("respIm", pg.math.median(respIm), min(respIm), max(respIm))
        JC = pg.matrix.CMatrix()
        self.createJacobian_(modelC, u, JC)
        for i in range(JC.rows()):
            #JC[i] *= 1.0/(modelC*modelC) * k[i]
            JC[i] /= (modelC * modelC) / k[i]
        self._J.mat(0).copy(pg.math.real(JC))
        self._J.mat(1).copy(pg.math.imag(JC))
        #self.createJacobian_(modelRe*0.0+1.0, pg.math.real(u), self._J.mat(1))
        #self.createJacobian_(modelRe*0.0+1.0, pg.math.imag(u), self._J.mat(2))
        #self.createJacobian_(modelRe*0.0+1.0, pg.math.imag(u), self._J.mat(3))
        sumsens0 = pg.Vector(self._J.mat(0).rows())
        sumsens1 = pg.Vector(self._J.mat(0).rows())
        sumsens2 = pg.Vector(self._J.mat(0).rows())
        for i in range(self._J.mat(0).rows()):
            #self._J.mat(0)[i] *= 1./modelRe / respRe[i]
            #self._J.mat(1)[i] *= 1./modelIm / respRe[i]
            #self._J.mat(2)[i] *= 1./modelRe / respIm[i]
            #self._J.mat(3)[i] *= 1./modelIm / respIm[i]
            #self._J.mat(0)[i] *= 1./(modelRe * modelRe) * k[i]
            #self._J.mat(1)[i] *= 1./(modelRe * modelIm) * k[i]
            #self._J.mat(2)[i] *= 1./(modelIm * modelRe) * k[i]
            #self._J.mat(3)[i] *= 1./(modelIm * modelIm) * k[i]
            sumsens0[i] = sum(self._J.mat(0)[i])
            sumsens1[i] = sum(self._J.mat(1)[i])
            sumsens2[i] = abs(sum(JC[i]))
        print(pg.math.median(sumsens0), min(sumsens0), max(sumsens0))
        print(pg.math.median(sumsens1), min(sumsens1), max(sumsens1))
        print(pg.math.median(sumsens2), min(sumsens2), max(sumsens2))
        self.data().set('k', k)
        self._J.recalcMatrixSize()
    else:
        # self.setVerbose(True)
        u = self.prepareJacobian_(model)
        #J = pg.Matrix()
        if self._J.rows() == 0:
            print('#' * 100)
            M1 = pg.Matrix()
            Jid = self._J.addMatrix(M1)
            self._J.addMatrixEntry(Jid, 0, 0)
        else:
            self._J.clean()
        self.createJacobian_(model, u, self._J.mat(0))
        self._J.recalcMatrixSize()