1# Copyright 2014-2020 The PySCF Developers. All Rights Reserved.
2#
3# Licensed under the Apache License, Version 2.0 (the "License");
4# you may not use this file except in compliance with the License.
5# You may obtain a copy of the License at
6#
7#     http://www.apache.org/licenses/LICENSE-2.0
8#
9# Unless required by applicable law or agreed to in writing, software
10# distributed under the License is distributed on an "AS IS" BASIS,
11# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12# See the License for the specific language governing permissions and
13# limitations under the License.
14
15'''
16Coupled Cluster
17===============
18
19Simple usage::
20
21    >>> from pyscf import gto, scf, cc
22    >>> mol = gto.M(atom='H 0 0 0; H 0 0 1')
23    >>> mf = scf.RHF(mol).run()
24    >>> cc.CCSD(mf).run()
25
26:func:`cc.CCSD` returns an instance of CCSD class.  Followings are parameters
27to control CCSD calculation.
28
29    verbose : int
30        Print level.  Default value equals to :class:`Mole.verbose`
31    max_memory : float or int
32        Allowed memory in MB.  Default value equals to :class:`Mole.max_memory`
33    conv_tol : float
34        converge threshold.  Default is 1e-7.
35    conv_tol_normt : float
36        converge threshold for norm(t1,t2).  Default is 1e-5.
37    max_cycle : int
38        max number of iterations.  Default is 50.
39    diis_space : int
40        DIIS space size.  Default is 6.
41    diis_start_cycle : int
42        The step to start DIIS.  Default is 0.
43    direct : bool
44        AO-direct CCSD. Default is False.
45    async_io : bool
46        Allow for asynchronous function execution. Default is True.
47    incore_complete : bool
48        Avoid all I/O. Default is False.
49    frozen : int or list
50        If integer is given, the inner-most orbitals are frozen from CC
51        amplitudes.  Given the orbital indices (0-based) in a list, both
52        occupied and virtual orbitals can be frozen in CC calculation.
53
54
55Saved results
56
57    converged : bool
58        CCSD converged or not
59    e_tot : float
60        Total CCSD energy (HF + correlation)
61    t1, t2 :
62        t1[i,a], t2[i,j,a,b]  (i,j in occ, a,b in virt)
63    l1, l2 :
64        Lambda amplitudes l1[i,a], l2[i,j,a,b]  (i,j in occ, a,b in virt)
65'''
66
67from pyscf.cc import ccsd
68from pyscf.cc import ccsd_lambda
69from pyscf.cc import ccsd_rdm
70from pyscf.cc import addons
71from pyscf.cc import rccsd
72from pyscf.cc import uccsd
73from pyscf.cc import gccsd
74from pyscf.cc import eom_rccsd
75from pyscf.cc import eom_uccsd
76from pyscf.cc import eom_gccsd
77from pyscf.cc import qcisd
78from pyscf import scf
79
80def CCSD(mf, frozen=None, mo_coeff=None, mo_occ=None):
81    if isinstance(mf, scf.uhf.UHF):
82        return UCCSD(mf, frozen, mo_coeff, mo_occ)
83    elif isinstance(mf, scf.ghf.GHF):
84        return GCCSD(mf, frozen, mo_coeff, mo_occ)
85    else:
86        return RCCSD(mf, frozen, mo_coeff, mo_occ)
87CCSD.__doc__ = ccsd.CCSD.__doc__
88
89scf.hf.SCF.CCSD = CCSD
90
91
92def RCCSD(mf, frozen=None, mo_coeff=None, mo_occ=None):
93    import numpy
94    from pyscf import lib
95    from pyscf.soscf import newton_ah
96    from pyscf.cc import dfccsd
97
98    if isinstance(mf, scf.uhf.UHF):
99        raise RuntimeError('RCCSD cannot be used with UHF method.')
100    elif isinstance(mf, scf.rohf.ROHF):
101        lib.logger.warn(mf, 'RCCSD method does not support ROHF method. ROHF object '
102                        'is converted to UHF object and UCCSD method is called.')
103        mf = scf.addons.convert_to_uhf(mf)
104        return UCCSD(mf, frozen, mo_coeff, mo_occ)
105
106    if isinstance(mf, newton_ah._CIAH_SOSCF) or not isinstance(mf, scf.hf.RHF):
107        mf = scf.addons.convert_to_rhf(mf)
108
109    if getattr(mf, 'with_df', None):
110        return dfccsd.RCCSD(mf, frozen, mo_coeff, mo_occ)
111
112    elif numpy.iscomplexobj(mo_coeff) or numpy.iscomplexobj(mf.mo_coeff):
113        return rccsd.RCCSD(mf, frozen, mo_coeff, mo_occ)
114
115    else:
116        return ccsd.CCSD(mf, frozen, mo_coeff, mo_occ)
117RCCSD.__doc__ = ccsd.CCSD.__doc__
118
119
120def UCCSD(mf, frozen=None, mo_coeff=None, mo_occ=None):
121    from pyscf.soscf import newton_ah
122
123    if isinstance(mf, newton_ah._CIAH_SOSCF) or not isinstance(mf, scf.uhf.UHF):
124        mf = scf.addons.convert_to_uhf(mf)
125
126    if getattr(mf, 'with_df', None):
127        raise NotImplementedError('DF-UCCSD')
128    else:
129        return uccsd.UCCSD(mf, frozen, mo_coeff, mo_occ)
130UCCSD.__doc__ = uccsd.UCCSD.__doc__
131
132
133def GCCSD(mf, frozen=None, mo_coeff=None, mo_occ=None):
134    from pyscf.soscf import newton_ah
135
136    if isinstance(mf, newton_ah._CIAH_SOSCF) or not isinstance(mf, scf.ghf.GHF):
137        mf = scf.addons.convert_to_ghf(mf)
138
139    if getattr(mf, 'with_df', None):
140        raise NotImplementedError('DF-GCCSD')
141    else:
142        return gccsd.GCCSD(mf, frozen, mo_coeff, mo_occ)
143GCCSD.__doc__ = gccsd.GCCSD.__doc__
144
145
146def QCISD(mf, frozen=None, mo_coeff=None, mo_occ=None):
147    if isinstance(mf, scf.uhf.UHF):
148        raise NotImplementedError
149    elif isinstance(mf, scf.ghf.GHF):
150        raise NotImplementedError
151    else:
152        return RQCISD(mf, frozen, mo_coeff, mo_occ)
153QCISD.__doc__ = qcisd.QCISD.__doc__
154
155scf.hf.SCF.QCISD = QCISD
156
157def RQCISD(mf, frozen=None, mo_coeff=None, mo_occ=None):
158    import numpy
159    from pyscf import lib
160    from pyscf.soscf import newton_ah
161
162    if isinstance(mf, scf.uhf.UHF):
163        raise RuntimeError('RQCISD cannot be used with UHF method.')
164    elif isinstance(mf, scf.rohf.ROHF):
165        lib.logger.warn(mf, 'RQCISD method does not support ROHF method. ROHF object '
166                        'is converted to UHF object and UQCISD method is called.')
167        mf = scf.addons.convert_to_uhf(mf)
168        raise NotImplementedError
169
170    if isinstance(mf, newton_ah._CIAH_SOSCF) or not isinstance(mf, scf.hf.RHF):
171        mf = scf.addons.convert_to_rhf(mf)
172
173    elif numpy.iscomplexobj(mo_coeff) or numpy.iscomplexobj(mf.mo_coeff):
174        raise NotImplementedError
175
176    else:
177        return qcisd.QCISD(mf, frozen, mo_coeff, mo_occ)
178RQCISD.__doc__ = qcisd.QCISD.__doc__
179
180
181def FNOCCSD(mf, thresh=1e-6, pct_occ=None, nvir_act=None):
182    """Frozen natural orbital CCSD
183
184    Attributes:
185        thresh : float
186            Threshold on NO occupation numbers. Default is 1e-6 (very conservative).
187        pct_occ : float
188            Percentage of total occupation number. Default is None. If present, overrides `thresh`.
189        nvir_act : int
190            Number of virtual NOs to keep. Default is None. If present, overrides `thresh` and `pct_occ`.
191    """
192    from pyscf import mp
193    pt = mp.MP2(mf).set(verbose=0).run()
194    frozen, no_coeff = pt.make_fno(thresh=thresh, pct_occ=pct_occ, nvir_act=nvir_act)
195    pt_no = mp.MP2(mf, frozen=frozen, mo_coeff=no_coeff).set(verbose=0).run()
196    mycc = ccsd.CCSD(mf, frozen=frozen, mo_coeff=no_coeff)
197    mycc.delta_emp2 = pt.e_corr - pt_no.e_corr
198    from pyscf.lib import logger
199    def _finalize(self):
200        '''Hook for dumping results and clearing up the object.'''
201        if self.converged:
202            logger.info(self, 'FNO-%s converged', self.__class__.__name__)
203        else:
204            logger.note(self, 'FNO-%s not converged', self.__class__.__name__)
205        logger.note(self, 'E(FNO-%s) = %.16g  E_corr = %.16g',
206                    self.__class__.__name__, self.e_tot, self.e_corr)
207        logger.note(self, 'E(FNO-%s+delta-MP2) = %.16g  E_corr = %.16g',
208                    self.__class__.__name__, self.e_tot+self.delta_emp2,
209                    self.e_corr+self.delta_emp2)
210        return self
211    mycc._finalize = _finalize.__get__(mycc, mycc.__class__)
212    return mycc
213