Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

importing gpu4pyscf breaks hessian #130

Open
Acetylsalicylsaeure opened this issue Mar 24, 2024 · 6 comments
Open

importing gpu4pyscf breaks hessian #130

Acetylsalicylsaeure opened this issue Mar 24, 2024 · 6 comments
Labels
bug Something isn't working

Comments

@Acetylsalicylsaeure
Copy link

i get an error when trying to calculate the hessian. it works fine with the usual code, but breaks when i import gpu4pyscf.

regular pySCF:

In [1]: from pyscf import gto
   ...: from pyscf.hessian import thermo
   ...: 
   ...: # First compute nuclear Hessian.
   ...: mol = gto.M(
   ...:     atom = '''O    0.   0.       0
   ...:               H    0.   -0.757   0.587
   ...:               H    0.    0.757   0.587''',
   ...:     basis = '631g')
   ...: 
   ...: mf = mol.RKS().run()
   ...: hessian = mf.Hessian().kernel()
converged SCF energy = -75.8179288833364

after importing gpu4pyscf:

In [1]: from gpu4pyscf.dft import rks

In [2]: from pyscf import gto
   ...: from pyscf.hessian import thermo
   ...: 
   ...: # First compute nuclear Hessian.
   ...: mol = gto.M(
   ...:     atom = '''O    0.   0.       0
   ...:               H    0.   -0.757   0.587
   ...:               H    0.    0.757   0.587''',
   ...:     basis = '631g')
   ...: 
   ...: mf = mol.RKS().run()
   ...: hessian = mf.Hessian().kernel()

converged SCF energy = -75.8179288833364
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In[2], line 12
      5 mol = gto.M(
      6     atom = '''O    0.   0.       0
      7               H    0.   -0.757   0.587
      8               H    0.    0.757   0.587''',
      9     basis = '631g')
     11 mf = mol.RKS().run()
---> 12 hessian = mf.Hessian().kernel()

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/hessian/rhf.py:524, in kernel(hessobj, mo_energy, mo_coeff, mo_occ, atmlst)
    521 else:
    522     hessobj.atmlst = atmlst
--> 524 de = hessobj.hess_elec(mo_energy, mo_coeff, mo_occ, atmlst=atmlst)
    525 hessobj.de = de.get() + hessobj.hess_nuc(hessobj.mol, atmlst=atmlst)
    526 mf = hessobj.base

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/hessian/rhf.py:59, in hess_elec(hessobj, mo_energy, mo_coeff, mo_occ, mo1, mo_e1, h1ao, atmlst, max_memory, verbose)
     57 mo_occ = cupy.asarray(mo_occ)
     58 mo_coeff = cupy.asarray(mo_coeff)
---> 59 de2 = hessobj.partial_hess_elec(mo_energy, mo_coeff, mo_occ, atmlst,
     60                                 max_memory, log)
     61 t1 = log.timer_debug1('hess elec', *t1)
     62 if h1ao is None:

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/hessian/rks.py:61, in partial_hess_elec(hessobj, mo_energy, mo_coeff, mo_occ, atmlst, max_memory, verbose)
     59 else:
     60     hyb = mf._numint.hybrid_coeff(mf.xc, spin=mol.spin)
---> 61 de2, ej, ek = rhf_hess._partial_hess_ejk(hessobj, mo_energy, mo_coeff, mo_occ,
     62                                          atmlst, max_memory, verbose,
     63                                          abs(hyb) > 1e-10)
     64 de2 += ej - hyb * ek  # (A,B,dR_A,dR_B)
     66 mem_now = lib.current_memory()[0]

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/hessian/rhf.py:181, in _partial_hess_ejk(hessobj, mo_energy, mo_coeff, mo_occ, atmlst, max_memory, verbose, with_k)
    178 vk1 = vk1.reshape(3,3,nao,nao)
    179 t1 = log.timer_debug1('contracting int2e_ipvip1 for atom %d'%ia, *t1)
--> 181 ej[i0,i0] += contract('xypq,pq->xy', vj1_diag[:,:,p0:p1], dm0[p0:p1])*2
    182 ek[i0,i0] += contract('xypq,pq->xy', vk1_diag[:,:,p0:p1], dm0[p0:p1])
    183 e1[i0,i0] -= contract('xypq,pq->xy', s1aa[:,:,p0:p1], dme0[p0:p1])*2

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/lib/cutensor.py:151, in contract(pattern, a, b, alpha, beta, out)
    146 def contract(pattern, a, b, alpha=1.0, beta=0.0, out=None):
    147     '''
    148     a wrapper for general tensor contraction
    149     pattern has to be a standard einsum notation
    150     '''
--> 151     return contraction(pattern, a, b, alpha, beta, out=out)

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/lib/cutensor.py:109, in contraction(pattern, a, b, alpha, beta, out, op_a, op_b, op_c, algo, jit_mode, compute_desc, ws_pref)
    105 alpha = np.asarray(alpha)
    106 beta = np.asarray(beta)
    108 cutensor_backend.contract(cutensor._get_handle().ptr, plan.ptr,
--> 109                          alpha.ctypes.data, a.data.ptr, b.data.ptr,
    110                          beta.ctypes.data, c.data.ptr, out.data.ptr,
    111                          ws.data.ptr, ws_size)
    113 return out

AttributeError: 'memoryview' object has no attribute 'ptr'

i get the same error after just calling to_gpu() method
@sunqm
Copy link
Collaborator

sunqm commented Mar 24, 2024

Please use gpu4pyscf.dft.RKS(mol) instead of mol.RKS(). The later creates a RKS object implemented in CPU. CPU module and GPU module are not compatible in this case.

There needs a check or a conversion on the ._scf attribute, when initializing hessian, gradients and other modules that rely on mf objects. @wxj6000

@wxj6000
Copy link
Collaborator

wxj6000 commented Mar 24, 2024

Please use gpu4pyscf.dft.RKS(mol) instead of mol.RKS(). The later creates a RKS object implemented in CPU. CPU module and GPU module are not compatible in this case.

There needs a check or a conversion on the ._scf attribute, when initializing hessian, gradients and other modules that rely on mf objects. @wxj6000

This is one of the cases I didn't check the API compatibility carefully. I will work on it this coming week.

@Acetylsalicylsaeure Thank you for your feedback. Please use the solution by Qiming at this point. We are in the process of improving the compatibility of PySCF and GPU4PySCF.

@Acetylsalicylsaeure
Copy link
Author

i was more worried about the occurence of pure pyscf code suddenly breaking due to gpu4pyscf being imported
using gpu4pyscf.dft.rks.RKS(mol) still leads to the same error

In [1]: import numpy as np
   ...: import pyscf
   ...: from pyscf import lib, df
   ...: from gpu4pyscf.dft import rks
   ...: lib.num_threads(8)
   ...: 
   ...: atom ='''
   ...: O       0.0000000000    -0.0000000000     0.1174000000
   ...: H      -0.7570000000    -0.0000000000    -0.4696000000
   ...: H       0.7570000000     0.0000000000    -0.4696000000
   ...: '''
   ...: 
   ...: mol = pyscf.M(atom=atom, basis='def2-tzvpp', verbose=0)
   ...: 

In [2]: mf = rks.RKS(mol)

In [3]: mf.kernel()
Out[3]: -75.90464116157013

In [4]: hessobj = mf.Hessian()

In [5]: hess = hessobj.kernel()
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In[5], line 1
----> 1 hess = hessobj.kernel()

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/hessian/rhf.py:524, in kernel(hessobj, mo_energy, mo_coeff, mo_occ, atmlst)
    521 else:
    522     hessobj.atmlst = atmlst
--> 524 de = hessobj.hess_elec(mo_energy, mo_coeff, mo_occ, atmlst=atmlst)
    525 hessobj.de = de.get() + hessobj.hess_nuc(hessobj.mol, atmlst=atmlst)
    526 mf = hessobj.base

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/hessian/rhf.py:59, in hess_elec(hessobj, mo_energy, mo_coeff, mo_occ, mo1, mo_e1, h1ao, atmlst, max_memory, verbose)
     57 mo_occ = cupy.asarray(mo_occ)
     58 mo_coeff = cupy.asarray(mo_coeff)
---> 59 de2 = hessobj.partial_hess_elec(mo_energy, mo_coeff, mo_occ, atmlst,
     60                                 max_memory, log)
     61 t1 = log.timer_debug1('hess elec', *t1)
     62 if h1ao is None:

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/hessian/rks.py:61, in partial_hess_elec(hessobj, mo_energy, mo_coeff, mo_occ, atmlst, max_memory, verbose)
     59 else:
     60     hyb = mf._numint.hybrid_coeff(mf.xc, spin=mol.spin)
---> 61 de2, ej, ek = rhf_hess._partial_hess_ejk(hessobj, mo_energy, mo_coeff, mo_occ,
     62                                          atmlst, max_memory, verbose,
     63                                          abs(hyb) > 1e-10)
     64 de2 += ej - hyb * ek  # (A,B,dR_A,dR_B)
     66 mem_now = lib.current_memory()[0]

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/hessian/rhf.py:181, in _partial_hess_ejk(hessobj, mo_energy, mo_coeff, mo_occ, atmlst, max_memory, verbose, with_k)
    178 vk1 = vk1.reshape(3,3,nao,nao)
    179 t1 = log.timer_debug1('contracting int2e_ipvip1 for atom %d'%ia, *t1)
--> 181 ej[i0,i0] += contract('xypq,pq->xy', vj1_diag[:,:,p0:p1], dm0[p0:p1])*2
    182 ek[i0,i0] += contract('xypq,pq->xy', vk1_diag[:,:,p0:p1], dm0[p0:p1])
    183 e1[i0,i0] -= contract('xypq,pq->xy', s1aa[:,:,p0:p1], dme0[p0:p1])*2

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/lib/cutensor.py:151, in contract(pattern, a, b, alpha, beta, out)
    146 def contract(pattern, a, b, alpha=1.0, beta=0.0, out=None):
    147     '''
    148     a wrapper for general tensor contraction
    149     pattern has to be a standard einsum notation
    150     '''
--> 151     return contraction(pattern, a, b, alpha, beta, out=out)

File ~/miniconda3/envs/openmm/lib/python3.11/site-packages/gpu4pyscf/lib/cutensor.py:109, in contraction(pattern, a, b, alpha, beta, out, op_a, op_b, op_c, algo, jit_mode, compute_desc, ws_pref)
    105 alpha = np.asarray(alpha)
    106 beta = np.asarray(beta)
    108 cutensor_backend.contract(cutensor._get_handle().ptr, plan.ptr,
--> 109                          alpha.ctypes.data, a.data.ptr, b.data.ptr,
    110                          beta.ctypes.data, c.data.ptr, out.data.ptr,
    111                          ws.data.ptr, ws_size)
    113 return out

AttributeError: 'memoryview' object has no attribute 'ptr'

In [6]:

however, calling .density_fit() on the RKS object fixes the error, is this expected behavior?

@Acetylsalicylsaeure
Copy link
Author

Please use gpu4pyscf.dft.RKS(mol) instead of mol.RKS(). The later creates a RKS object implemented in CPU. CPU module and GPU module are not compatible in this case.
There needs a check or a conversion on the ._scf attribute, when initializing hessian, gradients and other modules that rely on mf objects. @wxj6000

This is one of the cases I didn't check the API compatibility carefully. I will work on it this coming week.

@Acetylsalicylsaeure Thank you for your feedback. Please use the solution by Qiming at this point. We are in the process of improving the compatibility of PySCF and GPU4PySCF.

ah, i was too slow, thanks! i also experimented around a bit, works well with RHF

@wxj6000
Copy link
Collaborator

wxj6000 commented Mar 24, 2024

@Acetylsalicylsaeure Hessian calculation of regular SCF is not supported on GPU yet. It is only supported for density fitting scheme.

@wxj6000 wxj6000 added the bug Something isn't working label Apr 2, 2024
@wxj6000 wxj6000 mentioned this issue Apr 2, 2024
Merged
@wxj6000
Copy link
Collaborator

wxj6000 commented Apr 2, 2024

@Acetylsalicylsaeure This is a bug in GPU4PySCF RKS hessian module. PySCF hessian has been redefined. The issue has been fixed in the above PR. Thank you for your feedback.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
bug Something isn't working
Projects
None yet
Milestone
No milestone
Development

No branches or pull requests
3 participants
@sunqm @wxj6000 @Acetylsalicylsaeure