      SUBROUTINE PLEM(Z, W, IERR)
C
C     WEIERSTRASS P-FUNCTION IN THE LEMNISCATIC CASE
C     FOR COMPLEX ARGUMENT WITH UNIT PERIOD PARALLELOGRAM
C
      COMPLEX Z, Z1, Z4, Z6, W
      REAL ZR, ZI
      INTEGER IERR, M, N
C
C     REDUCTION TO FUNDAMENTAL PARALLELOGRAM
C
      ZR = REAL(Z) + 0.5E0
      ZI = AIMAG(Z) + 0.5E0
      M = INT(ZR)
      N = INT(ZI)
      IF (ZR.LT.0E0) M = M - 1
      IF (ZI.LT.0E0) N = N - 1
      Z1 = Z - FLOAT(M) - (0E0,1E0)*FLOAT(N)
C
C     IF Z1=0 THEN Z COINCIDES WITH A LATTICE POINT.
C     THE LATTICE POINTS ARE POLES FOR P.
C
      W = Z1*Z1
      ZR = ABS(REAL(W)) + ABS(AIMAG(W))
      IF (ZR.NE.0E0) GO TO 10
      IERR = 1
      RETURN
C
C     EVALUATION OF P(Z1)
C
   10 IERR = 0
      Z4 = W*W
      Z6 = Z4*W
      W = 1E0/W + 4E0*W*(3E0+Z4)/(1E0-Z4)**2 +
     * W*((((((((-7.233108E-11*Z4+1.714197273E-8)*Z4-2.5369036492E-7)*
     * Z4-7.98710206868E-6)*Z4+6.4850606909737E-4)*Z4+7.39624629362938E-
     * 3)*Z4+2.012382768497244E-2)*Z4+7.1177297543136598E-1)*
     * Z4-2.54636399353830738E0)/((((((((5.1161516E-10*Z4+6.61289408E-9)
     * *Z4+4.4618987048E-7)*Z4-8.42694918892E-6)*Z4+4.42886829095E-6)*
     * Z4-4.22629935217101E-3)*Z4+2.577496871700433E-2)*
     * Z4+4.2359940482277074E-1)*Z4+1E0)
      RETURN
      END
