X2C big 1-norm and reference wavefunction dependency


I’m running into some problems with x2c. The 1-norm is very big and the scf result also depends on the reference wavefunction. Each run time the numbers can be different as well. I’m using decontracted 6-31g* here, which is not designed for relativistic calculations, but SARC basis set (uncontracted) for U atom also has the same problem. I imagine this has to do with the linear dependency problem or/and unprotected inversion LAPACK, but how should I fix it? Thanks!

Here’s an example input:
molecule h2o {
H 1 R
H 1 R 2 A

R = 2.0
A = 104.5
units bohr

set {
basis 6-31G*-decon
basis_relativistic 6-31g*-decon
scf_type pk
set relativistic x2c
en1, wfn1 = energy(‘scf’, return_wfn=True)

set relativistic no
en2, wfn2 = energy(‘scf’, return_wfn=True, reference_wfn=wfn1)

set relativistic x2c
en3, wfn3 = energy(‘scf’, return_wfn=True, reference_wfn=wfn2)

en4, wfn4 = energy(‘scf’, return_wfn=True, reference_wfn=wfn3)

(en1) The 1-norm of |H_X2C - H_Dirac| is: 41862.990477155596
(en3) The 1-norm of |H_X2C - H_Dirac| is: 41863.183105109230
(en4) The 1-norm of |H_X2C - H_Dirac| is: 41868.473475888735

en1 = -84575.5564468204102013
en2 = -2400.4387878319948868
en3 = -84575.6301333945448278
en4 = -84575.4272162587294588

I believe any X2C questions go to @francesco.