Energy jump in relaxation of PbS supercell containing one vacancy

[murphy]

Dear all,

I am doing structural relaxation of a rock salt (but in close-packed form, hexagonal) 107-atom PbS supercell containing one Pb vacancy. The relaxation steps are as follows:
1. ISIF=2, KPOINTS=Gamma-point only
2. ISIF=2, KPOINTS=Gamma-centered 2x2x2
3. ISIF=2, KPOINTS=Gamma-centered 4x4x4
4. ISIF=3, KPOINTS=Gamman-centered 4x4x4

But the converged energy for each steps are jumping:
1. -233.27933023
2. -232.28252864
3. -232.07779428
4. -232.06740704

In another case, I found the energies are STEP 1< STEP 2< STEP 3< STEP4, which is really strange. Because how come the the energy with poorer settings is lower than the one from a higher-setting calculation?

Hopefully some experienced users can solve my puzzle. Thank you in advance!
PS. The INCAR for STEP 1 is attached below. Please inform me if you need more information.

Regards,
Murphy
=====================================================
INCAR for STEP 1:
# Startparameter for this run:
PREC = high
NPAR = 4

# Electronic Relaxation 1
ENCUT = 500.0
ENAUG = 700.0
EDIFF = 0.00001
LREAL = A
NELM = 60
NELMIN = 2
NELMDL = 6 # of ELM steps

# Ionic relaxation
EDIFFG = -0.01 #stopping-criterion for IOM
NSW = 60 #number of steps for IOM
IBRION = 2 #ionic relax: 0-MD 1-quasi-New 2-CG
ISIF = 2
ISYM = 1 #0-nonsym 1-usesym 2-fastsym
POTIM = 0.5000 #time-step for ionic-motion

#DOS related values:
ISMEAR = 0
SIGMA = 0.05 #broadening in eV -4-tet -1-fermi 0-gaus

# Electronic relaxation 2 (details)
IALGO = 48 #algorithm

LMAXMIX = 4

# Dimension of arrays
NBANDS = 720

[admin]

Energies from different settings are not comparable. Higher precision does not mean that the total energy is more negative.