docs: def init1(self, dae):
          """New initialization function"""
          self.servcall(dae)
          retval = True
  
          mva = self.system.mva
          self.p0 = mul(self.p0, self.gammap)
          self.q0 = mul(self.q0, self.gammaq)
  
          dae.y[self.vsd] = mul(dae.y[self.v], -sin(dae.y[self.a]))
          dae.y[self.vsq] = mul(dae.y[self.v], cos(dae.y[self.a]))
  
          rs = matrix(self.rs)
          rr = matrix(self.rr)
          xmu = matrix(self.xmu)
          x1 = matrix(self.xs) + xmu
          x2 = matrix(self.xr) + xmu
          Pg = matrix(self.p0)
          Qg = matrix(self.q0)
          Vc = dae.y[self.v]
          vsq = dae.y[self.vsq]
          vsd = dae.y[self.vsd]
  
          toSn = div(mva, self.Sn)  # to machine base
          toSb = self.Sn / mva  # to system base
  
          # rotor speed
          omega = 1 * (ageb(mva * Pg, self.Sn)) + \
              mul(0.5 + 0.5 * mul(Pg, toSn),
                  aandb(agtb(Pg, 0), altb(mva * Pg, self.Sn))) + \
              0.5 * (aleb(mva * Pg, 0))
  
          slip = 1 - omega
          theta = mul(self.Kp, mround(1000 * (omega - 1)) / 1000)
          theta = mmax(theta, 0)
  
          # prepare for the iterations
  
          irq = mul(-x1, toSb, (2 * omega - 1), div(1, Vc), div(1, xmu),
                    div(1, omega))
          isd = zeros(*irq.size)
          isq = zeros(*irq.size)
  
          # obtain ird isd isq
          for i in range(self.n):
              A = sparse([[-rs[i], vsq[i]], [x1[i], -vsd[i]]])
              B = matrix([vsd[i] - xmu[i] * irq[i], Qg[i]])
              linsolve(A, B)
              isd[i] = B[0]
              isq[i] = B[1]
          ird = -div(vsq + mul(rs, isq) + mul(x1, isd), xmu)
          vrd = -mul(rr, ird) + mul(
              slip,
              mul(x2, irq) + mul(xmu, isq))  # todo: check x1 or x2
          vrq = -mul(rr, irq) - mul(slip, mul(x2, ird) + mul(xmu, isd))
  
          # main iterations
          for i in range(self.n):
              mis = ones(6, 1)
              rows = [0, 0, 0, 1, 1, 1, 2, 2, 3, 3, 4, 4, 5]
              cols = [0, 1, 3, 0, 1, 2, 2, 4, 3, 5, 0, 1, 2]
  
              x = matrix([isd[i], isq[i], ird[i], irq[i], vrd[i], vrq[i]])
              # vals = [-rs, x1, xmu, -x1, -rs, -xmu, -rr,
              #         -1, -rr, -1, vsd, vsq, -xmu * Vc / x1]
  
              vals = [
                  -rs[i], x1[i], xmu[i], -x1[i], -rs[i], -xmu[i], -rr[i], -1,
                  -rr[i], -1, vsd[i], vsq[i], -xmu[i] * Vc[i] / x1[i]
              ]
              jac0 = spmatrix(vals, rows, cols, (6, 6), 'd')
              iter = 0
  
              while max(abs(mis)) > self.system.tds.config.tol:
                  if iter > 20:
                      logger.error(
                          'Initialization of DFIG <{}> failed.'.format(
                              self.name[i]))
                      retval = False
                      break
  
                  mis[0] = -rs[i] * x[0] + x1[i] * x[1] + xmu[i] * x[3] - vsd[i]
                  mis[1] = -rs[i] * x[1] - x1[i] * x[0] - xmu[i] * x[2] - vsq[i]
                  mis[2] = -rr[i] * x[2] + slip[i] * (
                      x2[i] * x[3] + xmu[i] * x[1]) - x[4]
                  mis[3] = -rr[i] * x[3] - slip[i] * (
                      x2[i] * x[2] + xmu[i] * x[0]) - x[5]
                  mis[4] = vsd[i] * x[0] + vsq[i] * x[1] + x[4] * x[2] + \
                      x[5] * x[3] - Pg[i]
                  mis[5] = -xmu[i] * Vc[i] * x[2] / x1[i] - \
                      Vc[i] * Vc[i] / x1[i] - Qg[i]
  
                  rows = [2, 2, 3, 3, 4, 4, 4, 4]
                  cols = [1, 3, 0, 2, 2, 3, 4, 5]
                  vals = [
                      slip[i] * xmu[i], slip[i] * x2[i], -slip[i] * xmu[i],
                      -slip[i] * x2[i], x[4], x[5], x[2], x[3]
                  ]
  
                  jac = jac0 + spmatrix(vals, rows, cols, (6, 6), 'd')
  
                  linsolve(jac, mis)
  
                  x -= mis
                  iter += 1
  
              isd[i] = x[0]
              isq[i] = x[1]
              ird[i] = x[2]
              irq[i] = x[3]
              vrd[i] = x[4]
              vrq[i] = x[5]
  
          dae.x[self.ird] = mul(self.u0, ird)
          dae.x[self.irq] = mul(self.u0, irq)
          dae.y[self.isd] = isd
          dae.y[self.isq] = isq
          dae.y[self.vrd] = vrd
          dae.y[self.vrq] = vrq
  
          dae.x[self.omega_m] = mul(self.u0, omega)
          dae.x[self.theta_p] = mul(self.u0, theta)
          dae.y[self.pwa] = mmax(mmin(2 * dae.x[self.omega_m] - 1, 1), 0)
  
          self.vref0 = mul(
              aneb(self.KV, 0), Vc - div(ird + div(Vc, xmu), self.KV))
          dae.y[self.vref] = self.vref0
          k = mul(div(x1, Vc, xmu, omega), toSb)
  
          self.irq_off = -mul(k, mmax(mmin(2 * omega - 1, 1), 0)) - irq
  
          # electrical torque in pu
          te = mul(
              xmu,
              mul(dae.x[self.irq], dae.y[self.isd]) - mul(
                  dae.x[self.ird], dae.y[self.isq]))
  
          for i in range(self.n):
              if te[i] < 0:
                  logger.error(
                      'Pe < 0 on bus <{}>. Wind speed initialize failed.'.
                      format(self.bus[i]))
                  retval = False
  
          # wind power in pu
          pw = mul(te, omega)
          dae.y[self.pw] = pw
  
          # wind speed initialization loop
  
          R = 4 * pi * self.system.freq * mul(self.R, self.ngb,
                                              div(1, self.npole))
          AA = pi * self.R**2
          vw = 0.9 * self.Vwn
  
          for i in range(self.n):
              mis = 1
              iter = 0
              while abs(mis) > self.system.tds.config.tol:
                  if iter > 50:
                      logger.error(
                          'Wind <{}> init failed. '
                          'Try increasing the nominal wind speed.'.
                          format(self.wind[i]))
                      retval = False
                      break
  
                  pw_iter, jac = self.windpower(self.ngen[i], self.rho[i], vw[i],
                                                AA[i], R[i], omega[i], theta[i])
  
                  mis = pw_iter - pw[i]
                  inc = -mis / jac[1]
                  vw[i] += inc
                  iter += 1
  
          # set wind speed
          dae.x[self.vw] = div(vw, self.Vwn)
  
          lamb = div(omega, vw, div(1, R))
          ilamb = div(1,
                      (div(1, lamb + 0.08 * theta) - div(0.035, theta**3 + 1)))
          cp = 0.22 * mul(
              div(116, ilamb) - 0.4 * theta - 5, exp(div(-12.5, ilamb)))
  
          dae.y[self.lamb] = lamb
          dae.y[self.ilamb] = ilamb
          dae.y[self.cp] = cp
  
          self.system.rmgen(self.gen)
  
          if not retval:
              logger.error('DFIG initialization failed')
  
          return retval