adc-connect/adcc | Build 24657644018 | adcc/adc_pp/state_diffdm_2p.py | Coveralls

1

#!/usr/bin/env python3

2

## vi: tabstop=4 shiftwidth=4 softtabstop=4 expandtab

3

## ---------------------------------------------------------------------

4

##

5

## Copyright (C) 2026 by the adcc authors

6

##

7

## This file is part of adcc.

8

##

9

## adcc is free software: you can redistribute it and/or modify

10

## it under the terms of the GNU General Public License as published

11

## by the Free Software Foundation, either version 3 of the License, or

12

## (at your option) any later version.

13

##

14

## adcc is distributed in the hope that it will be useful,

15

## but WITHOUT ANY WARRANTY; without even the implied warranty of

16

## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the

17

## GNU General Public License for more details.

18

##

19

## You should have received a copy of the GNU General Public License

20

## along with adcc. If not, see <http://www.gnu.org/licenses/>.

21

##

22

## ---------------------------------------------------------------------

31

33

34

38

42

44

45

    # one-particle ISR(0) diffdm

48

50

        # N^4: O^4 / N^4: O^4

51

        - 4.0 * einsum("il,jk->ijkl", p1_oo, d_oo)

52

    ).antisymmetrise(0, 1).antisymmetrise(2, 3)

54

        # N^4: O^2V^2 / N^4: O^2V^2

55

        - 1.0 * einsum("ja,ib->iajb", u1, u1)

56

        # N^4: O^2V^2 / N^4: O^2V^2

57

        + 1.0 * einsum("ab,ij->iajb", p1_vv, d_oo)

58

60

61

65

69

71

    # TODO move to intermediates!

72

    # one-particle ISR(0) diffdm

74

75

    # ISR(2) ISR intermediate (TODO Move to intermediates)

77

    # new ones

79

81

        # N^4: O^2V^2 / N^4: O^2V^2

82

        + 4.0 * (

83

            einsum("ib,ja->ijab", ru1, u1)

84

        ).antisymmetrise(0, 1).antisymmetrise(2, 3)

85

        # N^5: O^3V^2 / N^4: O^2V^2

86

        + 2.0 * (

87

            einsum("ijka,kb->ijab", ru1_ooov, u1)

88

        ).antisymmetrise(2, 3)

89

        # N^5: O^3V^2 / N^4: O^2V^2

90

        - 2.0 * (

91

            einsum("jk,ikab->ijab", p1_oo, t2)

92

        ).antisymmetrise(0, 1)

93

94

96

97

101

103

        # ISR(1)-d

106

108

            # N^5: O^3V^2 / N^4: O^2V^2

109

            - 2.0 * einsum("kb,ijab->ijka", u1, u2)

110

111

113

            # N^5: O^2V^3 / N^4: O^1V^3

114

            - 2.0 * einsum("ja,ijbc->iabc", u1, u2)

115

117

        # no doubles contribution

120

121

129

133

135

        # N^6: O^4V^2 / N^4: O^2V^2

136

        + 2.0 * einsum("ijab,klab->ijkl", u2, u2)

137

        # N^6: O^5V^1 / N^4: O^2V^2

138

        - 1.0 * einsum("klmc,ijmc->ijkl", einsum("mb,klbc->klmc", u1, t2),

139

                       einsum("ma,ijac->ijmc", u1, t2))

140

+ (

141

            # N^5: O^3V^2 / N^4: O^2V^2

142

            - 8.0 * einsum("ik,jl->ijkl", einsum("imab,kmab->ik", u2, u2), d_oo)

143

            # N^4: O^4 / N^4: O^4

144

            + 4.0 * einsum("il,jk->ijkl", einsum("ia,la->il", u1, u1), p0.oo)

145

            # N^6: O^5V^1 / N^4: O^2V^2

146

            + 4.0 * einsum("jlmc,ikmc->ijkl", einsum("lb,jmbc->jlmc", u1, t2),

147

                           einsum("ia,kmac->ikmc", u1, t2))

148

            # N^4: O^2V^2 / N^4: O^2V^2

149

            + 4.0 * einsum("il,jk->ijkl",

150

                           einsum("lc,ic->il", einsum("nb,lnbc->lc", u1, t2),

151

                                  einsum("ma,imac->ic", u1, t2)), d_oo)

152

            # N^5: O^3V^2 / N^4: O^2V^2

153

            + 4.0 * einsum("ik,jl->ijkl", einsum("kmnc,imnc->ik",

154

                                                 einsum("nb,kmbc->kmnc", u1, t2),

155

                                                 einsum("na,imac->imnc", u1, t2)),

156

                           d_oo)

157

            # N^5: O^3V^2 / N^4: O^2V^2

158

            + 2.0 * einsum("ik,jl->ijkl",

159

                           einsum("inbc,knbc->ik",

160

                                  einsum("mn,imbc->inbc",

161

                                         einsum("ia,la->il", u1, u1), t2), t2),

162

                           d_oo)

163

        ).antisymmetrise(0, 1).antisymmetrise(2, 3)

164

+ (

165

            # N^5: O^3V^2 / N^4: O^2V^2

166

            + 4.0 * einsum("iklb,jb->ijkl",

167

                           einsum("ic,klbc->iklb", einsum("ma,imac->ic", u1, t2),

168

                                  t2), u1)

169

            # N^6: O^4V^2 / N^4: O^2V^2

170

            + 2.0 * einsum("jklm,im->ijkl",

171

                           einsum("jmbc,klbc->jklm", t2, t2),

172

                           einsum("ia,la->il", u1, u1))

173

        ).antisymmetrise(0, 1).symmetrise([(0, 2), (1, 3)])

174

+ (

175

            # N^4: O^4 / N^4: O^4

176

            + 4.0 * einsum("il,jk->ijkl", einsum("im,lm->il",

177

                                                 einsum("ia,la->il", u1, u1),

178

                                                 p0.oo), d_oo)

179

            # N^4: O^2V^2 / N^4: O^2V^2

180

            + 4.0 * einsum("ik,jl->ijkl",

181

                           einsum("kb,ib->ik",

182

                                  einsum("mc,kmbc->kb",

183

                                         einsum("na,mnac->mc", u1, t2), t2), u1),

184

                           d_oo)

185

        ).antisymmetrise(0, 1).antisymmetrise(2, 3).symmetrise([(0, 2), (1, 3)])

186

188

        # N^6: O^4V^2 / N^4: O^2V^2

189

        + 1.0 * einsum("ijkl,la->ijka", einsum("klbc,ijbc->ijkl", u2, t2), u1)

190

        # N^5: O^3V^2 / N^4: O^2V^2

191

        + 2.0 * einsum("kb,ijab->ijka", einsum("lc,klbc->kb", u1, u2), t2)

192

+ (

193

            # N^4: O^2V^2 / N^4: O^2V^2

194

            - 4.0 * einsum("ja,ik->ijka", einsum("lb,jlab->ja", u1, u2), d_oo)

195

            # N^5: O^3V^2 / N^4: O^2V^2

196

            + 2.0 * einsum("jk,ia->ijka", einsum("klbc,jlbc->jk", u2, t2), u1)

197

            # N^6: O^4V^2 / N^4: O^2V^2

198

            - 4.0 * einsum("jklb,ilab->ijka", einsum("jc,klbc->jklb", u1, u2), t2)

199

            # N^5: O^3V^2 / N^4: O^2V^2

200

            + 2.0 * einsum("ja,ik->ijka",

201

                           einsum("jlmb,lmab->ja",

202

                                  einsum("jc,lmbc->jlmb", u1, u2), t2), d_oo)

203

            # N^4: O^2V^2 / N^4: O^2V^2

204

            - 4.0 * einsum("ia,jk->ijka",

205

                           einsum("lb,ilab->ia",

206

                                  einsum("mc,lmbc->lb", u1, u2), t2), d_oo)

207

            # N^5: O^3V^2 / N^4: O^2V^2

208

            + 2.0 * einsum("ia,jk->ijka",

209

                           einsum("il,la->ia",

210

                                  einsum("lmbc,imbc->il", u2, t2), u1), d_oo)

211

            # N^4: O^3V^1 / N^4: O^3V^1

212

            + 2.0 * einsum("jk,ia->ijka", einsum("kb,jb->jk", u1, p0.ov), u1)

213

            # N^4: O^3V^1 / N^4: O^3V^1

214

            + 2.0 * einsum("jk,ia->ijka", einsum("jb,kb->jk", u1, u1), p0.ov)

215

            # N^4: O^3V^1 / N^4: O^3V^1

216

            + 2.0 * einsum("ia,jk->ijka",

217

                           einsum("il,la->ia",

218

                                  einsum("lb,ib->il", u1, p0.ov), u1), d_oo)

219

            # N^4: O^3V^1 / N^4: O^3V^1

220

            + 2.0 * einsum("ia,jk->ijka",

221

                           einsum("il,la->ia", einsum("ia,la->il", u1, u1),

222

                                  p0.ov), d_oo)

223

        ).antisymmetrise(0, 1)

224

226

        # N^4: O^2V^2 / N^4: O^2V^2

227

        + 4.0 * einsum("ib,ja->ijab",

228

                       einsum("kc,ikbc->ib", u1, td2),

229

                       u1).antisymmetrise(0, 1).antisymmetrise(2, 3)

230

        # N^5: O^3V^2 / N^4: O^2V^2

231

        + 2.0 * einsum("ijka,kb->ijab",

232

                       einsum("kc,ijac->ijka", u1, td2), u1).antisymmetrise(2, 3)

233

        # N^5: O^3V^2 / N^4: O^2V^2

234

        + 2.0 * einsum("jk,ikab->ijab", einsum("jc,kc->jk", u1, u1), td2)

235

237

        # N^6: O^3V^3 / N^4: O^2V^2

238

        - 4.0 * einsum("jkac,ikbc->iajb", u2, u2)

239

        # N^5: O^2V^3 / N^4: O^2V^2

240

        + 2.0 * einsum("ab,ij->iajb", einsum("klac,klbc->ab", u2, u2), d_oo)

241

        # N^6: O^4V^2 / N^4: O^2V^2

242

        + 1.0 * einsum("iklb,jkla->iajb", einsum("ld,ikbd->iklb", u1, t2),

243

                       einsum("lc,jkac->jkla", u1, t2))

244

        # N^4: O^2V^2 / N^4: O^2V^2

245

        + 1.0 * einsum("ab,ij->iajb", einsum("ka,kb->ab", u1, u1), p0.oo)

246

        # N^4: O^2V^2 / N^4: O^2V^2

247

        - 1.0 * einsum("ij,ab->iajb", einsum("ic,jc->ij", u1, u1), p0.vv)

248

        # N^6: O^3V^3 / N^4: O^2V^2

249

        + 1.0 * einsum("jlac,ilbc->iajb",

250

                       einsum("kl,jkac->jlac", einsum("kd,ld->kl", u1, u1), t2), t2)

251

        # N^6: O^4V^2 / N^4: O^2V^2

252

        + 0.5 * einsum("ijla,lb->iajb",

253

                       einsum("ijkl,ka->ijla",

254

                              einsum("ikcd,jlcd->ijkl", t2, t2), u1), u1)

255

        # N^6: O^4V^2 / N^4: O^2V^2

256

        + 0.5 * einsum("jklb,ikla->iajb",

257

                       einsum("jd,klbd->jklb", u1, t2),

258

                       einsum("ic,klac->ikla", u1, t2))

259

        # N^4: O^2V^2 / N^4: O^2V^2

260

        - 1.0 * einsum("ib,ja->iajb",

261

                       einsum("kd,ikbd->ib", u1, t2), einsum("lc,jlac->ja", u1, t2))

262

        # N^5: O^2V^3 / N^4: O^2V^2

263

        + 1.0 * einsum("ab,ij->iajb",

264

                       einsum("lmad,lmbd->ab",

265

                              einsum("km,klad->lmad",

266

                                     einsum("kc,mc->km", u1, u1), t2), t2), d_oo)

267

        # N^4: O^2V^2 / N^4: O^2V^2

268

        - 1.0 * einsum("ab,ij->iajb",

269

                       einsum("lb,la->ab",

270

                              einsum("kd,klbd->lb", u1, t2),

271

                              einsum("mc,lmac->la", u1, t2)), d_oo)

272

        # N^5: O^3V^2 / N^4: O^2V^2

273

        - 0.5 * einsum("ab,ij->iajb",

274

                       einsum("klmb,klma->ab",

275

                              einsum("md,klbd->klmb", u1, t2),

276

                              einsum("mc,klac->klma", u1, t2)), d_oo)

277

+ (

278

            # N^4: O^2V^2 / N^4: O^2V^2

279

            + 1.0 * einsum("ib,ja->iajb", einsum("ic,bc->ib", u1, p0.vv), u1)

280

            # N^4: O^2V^2 / N^4: O^2V^2

281

            - 1.0 * einsum("ib,ja->iajb", einsum("kb,ik->ib", u1, p0.oo), u1)

282

            # N^5: O^3V^2 / N^4: O^2V^2

283

            + 2.0 * einsum("ijlb,la->iajb",

284

                           einsum("jc,ilbc->ijlb",

285

                                  einsum("kd,jkcd->jc", u1, t2), t2), u1)

286

            # N^5: O^3V^2 / N^4: O^2V^2

287

            + 2.0 * einsum("kb,ijka->iajb",

288

                           einsum("ld,klbd->kb", u1, t2),

289

                           einsum("ic,jkac->ijka", u1, t2))

290

            # N^6: O^4V^2 / N^4: O^2V^2

291

            - 2.0 * einsum("ijkb,ka->iajb",

292

                           einsum("jklc,ilbc->ijkb",

293

                                  einsum("kd,jlcd->jklc", u1, t2), t2), u1)

294

            # N^6: O^4V^2 / N^4: O^2V^2

295

            - 2.0 * einsum("ijlb,la->iajb",

296

                           einsum("ijkc,klbc->ijlb",

297

                                  einsum("id,jkcd->ijkc", u1, t2), t2), u1)

298

            # N^6: O^3V^3 / N^4: O^2V^2

299

            - 2.0 * einsum("ikbc,jkac->iajb",

300

                           einsum("il,klbc->ikbc",

301

                                  einsum("id,ld->il", u1, u1), t2), t2)

302

            # N^4: O^2V^2 / N^4: O^2V^2

303

            - 1.0 * einsum("ib,ja->iajb",

304

                           einsum("kc,ikbc->ib",

305

                                  einsum("ld,klcd->kc", u1, t2), t2), u1)

306

            # N^4: O^2V^2 / N^4: O^2V^2

307

            - 1.0 * einsum("ab,ij->iajb",

308

                           einsum("kb,ka->ab",

309

                                  einsum("kc,bc->kb", u1, p0.vv), u1), d_oo)

310

            # N^4: O^2V^2 / N^4: O^2V^2

311

            - 1.0 * einsum("ab,ij->iajb",

312

                           einsum("mb,ma->ab",

313

                                  einsum("ld,lmbd->mb",

314

                                         einsum("kc,klcd->ld", u1, t2),

315

                                         t2), u1), d_oo)

316

        ).symmetrise([(0, 2), (1, 3)])

317

319

        # N^6: O^3V^3 / N^4: O^1V^3

320

        + 1.0 * einsum("ijka,jkbc->iabc", einsum("id,jkad->ijka", u1, u2), t2)

321

        # N^5: O^2V^3 / N^4: O^1V^3

322

        + 2.0 * einsum("ja,ijbc->iabc", einsum("kd,jkad->ja", u1, u2), t2)

323

+ (

324

            # N^5: O^2V^3 / N^4: O^1V^3

325

            + 2.0 * einsum("ac,ib->iabc", einsum("jkad,jkcd->ac", u2, t2), u1)

326

            # N^6: O^3V^3 / N^4: O^1V^3

327

            - 4.0 * einsum("ikab,kc->iabc", einsum("jkad,ijbd->ikab", u2, t2), u1)

328

            # N^4: O^1V^3 / N^4: O^1V^3

329

            + 2.0 * einsum("ac,ib->iabc", einsum("ja,jc->ac", u1, p0.ov), u1)

330

            # N^4: O^1V^3 / N^4: O^1V^3

331

            + 2.0 * einsum("ac,ib->iabc", einsum("ja,jc->ac", u1, u1), p0.ov)

332

        ).antisymmetrise(2, 3)

333

335

        # N^6: O^2V^4 / N^4: V^4

336

        + 2 * einsum("ijab,ijcd->abcd", u2, u2)

337

        # N^6: O^2V^4 / N^4: V^4

338

        - 1 * einsum("jkab,jkcd->abcd",

339

                     einsum("ij,ikab->jkab", einsum("ie,je->ij", u1, u1), t2), t2)

340

+ (

341

            # N^4: V^4 / N^4: V^4

342

            + 4.0 * einsum("ac,bd->abcd", einsum("ia,ic->ac", u1, u1), p0.vv)

343

            # N^6: O^3V^3 / N^4: V^4

344

            + 4.0 * einsum("jabc,jd->abcd",

345

                           einsum("ijbc,ia->jabc",

346

                                  einsum("jkbe,ikce->ijbc", t2, t2), u1), u1)

347

        ).antisymmetrise(0, 1).antisymmetrise(2, 3)

348

+ (

349

            # N^5: O^1V^4 / N^4: V^4

350

            + 4.0 * einsum("ka,kbcd->abcd",

351

                           einsum("ie,ikae->ka", u1, t2),

352

                           einsum("jb,jkcd->kbcd", u1, t2))

353

            # N^6: O^3V^3 / N^4: V^4

354

            + 2.0 * einsum("ibcd,ia->abcd",

355

                           einsum("ijkb,jkcd->ibcd",

356

                                  einsum("ie,jkbe->ijkb", u1, t2), t2), u1)

357

        ).antisymmetrise(0, 1).symmetrise([(0, 2), (1, 3)])

358

360

361

362

# dict controlling the dispatch of the state_diffdm function

364

    "isr0": diffdm_isr0_2p,

365

    "isr1s": diffdm_isr1s_2p,

366

    "isr1": diffdm_isr1_2p,

367

    "isr2": diffdm_isr2_2p,

368

369

370

372

"""

373

    Compute the two-particle difference density matrix of an excited state

374

    in the MO basis.

375

376

    Parameters

377

    ----------

378

    method : str, IsrMethod

379

        The method to use for the computation (e.g. "isr2")

380

    ground_state : LazyMp

381

        The ground state upon which the excitation was based

382

    amplitude : AmplitudeVector

383

        The amplitude vector

384

    intermediates : adcc.Intermediates

385

        Intermediates from the ADC calculation to reuse

386

"""

395

398

                                  f"for {method.name}.")

399

    else:

adc-connect / adcc / 24657644018

Source File
Press 'n' to go to next uncovered line, 'b' for previous

1	#!/usr/bin/env python3
2	## vi: tabstop=4 shiftwidth=4 softtabstop=4 expandtab
3	## ---------------------------------------------------------------------
4	##
5	## Copyright (C) 2026 by the adcc authors
6	##
7	## This file is part of adcc.
8	##
9	## adcc is free software: you can redistribute it and/or modify
10	## it under the terms of the GNU General Public License as published
11	## by the Free Software Foundation, either version 3 of the License, or
12	## (at your option) any later version.
13	##
14	## adcc is distributed in the hope that it will be useful,
15	## but WITHOUT ANY WARRANTY; without even the implied warranty of
16	## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17	## GNU General Public License for more details.
18	##
19	## You should have received a copy of the GNU General Public License
20	## along with adcc. If not, see <http://www.gnu.org/licenses/>.
21	##
22	## ---------------------------------------------------------------------
23	from adcc import block as b	2✔
24	from adcc.LazyMp import LazyMp	2✔
25	from adcc.AdcMethod import IsrMethod	2✔
26	from adcc.functions import einsum, zeros_like	2✔
27	from adcc.Intermediates import Intermediates	2✔
28	from adcc.AmplitudeVector import AmplitudeVector	2✔
29	from adcc.TwoParticleDensity import TwoParticleDensity	2✔
30	from adcc.NParticleOperator import OperatorSymmetry	2✔
31
32	from .util import check_doubles_amplitudes, check_singles_amplitudes	2✔
33
34
35	def diffdm_isr0_2p(mp, amplitude, intermediates):	2✔
36	check_singles_amplitudes([b.o, b.v], amplitude)	2✔
37	u1 = amplitude.ph	2✔
38
39	hf = mp.reference_state	2✔
40	d_oo = zeros_like(hf.foo)	2✔
41	d_oo.set_mask("ii", 1)	2✔
42
43	dm = TwoParticleDensity(mp, symmetry=OperatorSymmetry.HERMITIAN)	2✔
44
45	# one-particle ISR(0) diffdm
46	p1_oo = -einsum("ia,la->il", u1, u1).evaluate()	2✔
47	p1_vv = einsum("ka,kb->ab", u1, u1).evaluate()	2✔
48
49	dm.oooo = (	2✔
50	# N^4: O^4 / N^4: O^4
51	- 4.0 * einsum("il,jk->ijkl", p1_oo, d_oo)
52	).antisymmetrise(0, 1).antisymmetrise(2, 3)
53	dm.ovov = (	2✔
54	# N^4: O^2V^2 / N^4: O^2V^2
55	- 1.0 * einsum("ja,ib->iajb", u1, u1)
56	# N^4: O^2V^2 / N^4: O^2V^2
57	+ 1.0 * einsum("ab,ij->iajb", p1_vv, d_oo)
58	)
59	return dm	2✔
60
61
62	def diffdm_isr1s_2p(mp, amplitude, intermediates):	2✔
63	dm = diffdm_isr0_2p(mp, amplitude, intermediates) # Get ISR(0) result	2✔
64	u1 = amplitude.ph	2✔
65
66	hf = mp.reference_state	2✔
67	d_oo = zeros_like(hf.foo)	2✔
68	d_oo.set_mask("ii", 1)	2✔
69
70	t2 = mp.t2(b.oovv)	2✔
71	# TODO move to intermediates!
72	# one-particle ISR(0) diffdm
73	p1_oo = -einsum("ia,la->il", u1, u1).evaluate()	2✔
74
75	# ISR(2) ISR intermediate (TODO Move to intermediates)
76	ru1 = einsum("ijab,jb->ia", t2, u1).evaluate()	2✔
77	# new ones
78	ru1_ooov = einsum("kc,ijac->ijka", u1, t2).evaluate()	2✔
79
80	dm.oovv += (	2✔
81	# N^4: O^2V^2 / N^4: O^2V^2
82	+ 4.0 * (
83	einsum("ib,ja->ijab", ru1, u1)
84	).antisymmetrise(0, 1).antisymmetrise(2, 3)
85	# N^5: O^3V^2 / N^4: O^2V^2
86	+ 2.0 * (
87	einsum("ijka,kb->ijab", ru1_ooov, u1)
88	).antisymmetrise(2, 3)
89	# N^5: O^3V^2 / N^4: O^2V^2
90	- 2.0 * (
91	einsum("jk,ikab->ijab", p1_oo, t2)
92	).antisymmetrise(0, 1)
93	)
94
95	return dm	2✔
96
97
98	def diffdm_isr1_2p(mp, amplitude, intermediates):	2✔
99	dm = diffdm_isr1s_2p(mp, amplitude, intermediates) # Get ISR(1)-s result	2✔
100	u1 = amplitude.ph	2✔
101
102	try:	2✔
103	# ISR(1)-d
104	check_doubles_amplitudes([b.o, b.o, b.v, b.v], amplitude)	2✔
105	u2 = amplitude.pphh	2✔
106
107	dm.ooov += (	2✔
108	# N^5: O^3V^2 / N^4: O^2V^2
109	- 2.0 * einsum("kb,ijab->ijka", u1, u2)
110	)
111
112	dm.ovvv += (	2✔
113	# N^5: O^2V^3 / N^4: O^1V^3
114	- 2.0 * einsum("ja,ijbc->iabc", u1, u2)
115	)
NEW 116	except ValueError:	×
117	# no doubles contribution
NEW 118	pass	×
119	return dm	2✔
120
121
122	def diffdm_isr2_2p(mp, amplitude, intermediates):	2✔
123	dm = diffdm_isr1_2p(mp, amplitude, intermediates) # Get ISR(1) result	2✔
124	check_doubles_amplitudes([b.o, b.o, b.v, b.v], amplitude)	2✔
125	u1, u2 = amplitude.ph, amplitude.pphh	2✔
126	hf = mp.reference_state	2✔
127	d_oo = zeros_like(hf.foo)	2✔
128	d_oo.set_mask("ii", 1)	2✔
129
130	t2 = mp.t2(b.oovv)	2✔
131	td2 = mp.td2(b.oovv)	2✔
132	p0 = mp.mp2_diffdm	2✔
133
134	dm.oooo += (	2✔
135	# N^6: O^4V^2 / N^4: O^2V^2
136	+ 2.0 * einsum("ijab,klab->ijkl", u2, u2)
137	# N^6: O^5V^1 / N^4: O^2V^2
138	- 1.0 * einsum("klmc,ijmc->ijkl", einsum("mb,klbc->klmc", u1, t2),
139	einsum("ma,ijac->ijmc", u1, t2))
140	+ (
141	# N^5: O^3V^2 / N^4: O^2V^2
142	- 8.0 * einsum("ik,jl->ijkl", einsum("imab,kmab->ik", u2, u2), d_oo)
143	# N^4: O^4 / N^4: O^4
144	+ 4.0 * einsum("il,jk->ijkl", einsum("ia,la->il", u1, u1), p0.oo)
145	# N^6: O^5V^1 / N^4: O^2V^2
146	+ 4.0 * einsum("jlmc,ikmc->ijkl", einsum("lb,jmbc->jlmc", u1, t2),
147	einsum("ia,kmac->ikmc", u1, t2))
148	# N^4: O^2V^2 / N^4: O^2V^2
149	+ 4.0 * einsum("il,jk->ijkl",
150	einsum("lc,ic->il", einsum("nb,lnbc->lc", u1, t2),
151	einsum("ma,imac->ic", u1, t2)), d_oo)
152	# N^5: O^3V^2 / N^4: O^2V^2
153	+ 4.0 * einsum("ik,jl->ijkl", einsum("kmnc,imnc->ik",
154	einsum("nb,kmbc->kmnc", u1, t2),
155	einsum("na,imac->imnc", u1, t2)),
156	d_oo)
157	# N^5: O^3V^2 / N^4: O^2V^2
158	+ 2.0 * einsum("ik,jl->ijkl",
159	einsum("inbc,knbc->ik",
160	einsum("mn,imbc->inbc",
161	einsum("ia,la->il", u1, u1), t2), t2),
162	d_oo)
163	).antisymmetrise(0, 1).antisymmetrise(2, 3)
164	+ (
165	# N^5: O^3V^2 / N^4: O^2V^2
166	+ 4.0 * einsum("iklb,jb->ijkl",
167	einsum("ic,klbc->iklb", einsum("ma,imac->ic", u1, t2),
168	t2), u1)
169	# N^6: O^4V^2 / N^4: O^2V^2
170	+ 2.0 * einsum("jklm,im->ijkl",
171	einsum("jmbc,klbc->jklm", t2, t2),
172	einsum("ia,la->il", u1, u1))
173	).antisymmetrise(0, 1).symmetrise([(0, 2), (1, 3)])
174	+ (
175	# N^4: O^4 / N^4: O^4
176	+ 4.0 * einsum("il,jk->ijkl", einsum("im,lm->il",
177	einsum("ia,la->il", u1, u1),
178	p0.oo), d_oo)
179	# N^4: O^2V^2 / N^4: O^2V^2
180	+ 4.0 * einsum("ik,jl->ijkl",
181	einsum("kb,ib->ik",
182	einsum("mc,kmbc->kb",
183	einsum("na,mnac->mc", u1, t2), t2), u1),
184	d_oo)
185	).antisymmetrise(0, 1).antisymmetrise(2, 3).symmetrise([(0, 2), (1, 3)])
186	)
187	dm.ooov += (	2✔
188	# N^6: O^4V^2 / N^4: O^2V^2
189	+ 1.0 * einsum("ijkl,la->ijka", einsum("klbc,ijbc->ijkl", u2, t2), u1)
190	# N^5: O^3V^2 / N^4: O^2V^2
191	+ 2.0 * einsum("kb,ijab->ijka", einsum("lc,klbc->kb", u1, u2), t2)
192	+ (
193	# N^4: O^2V^2 / N^4: O^2V^2
194	- 4.0 * einsum("ja,ik->ijka", einsum("lb,jlab->ja", u1, u2), d_oo)
195	# N^5: O^3V^2 / N^4: O^2V^2
196	+ 2.0 * einsum("jk,ia->ijka", einsum("klbc,jlbc->jk", u2, t2), u1)
197	# N^6: O^4V^2 / N^4: O^2V^2
198	- 4.0 * einsum("jklb,ilab->ijka", einsum("jc,klbc->jklb", u1, u2), t2)
199	# N^5: O^3V^2 / N^4: O^2V^2
200	+ 2.0 * einsum("ja,ik->ijka",
201	einsum("jlmb,lmab->ja",
202	einsum("jc,lmbc->jlmb", u1, u2), t2), d_oo)
203	# N^4: O^2V^2 / N^4: O^2V^2
204	- 4.0 * einsum("ia,jk->ijka",
205	einsum("lb,ilab->ia",
206	einsum("mc,lmbc->lb", u1, u2), t2), d_oo)
207	# N^5: O^3V^2 / N^4: O^2V^2
208	+ 2.0 * einsum("ia,jk->ijka",
209	einsum("il,la->ia",
210	einsum("lmbc,imbc->il", u2, t2), u1), d_oo)
211	# N^4: O^3V^1 / N^4: O^3V^1
212	+ 2.0 * einsum("jk,ia->ijka", einsum("kb,jb->jk", u1, p0.ov), u1)
213	# N^4: O^3V^1 / N^4: O^3V^1
214	+ 2.0 * einsum("jk,ia->ijka", einsum("jb,kb->jk", u1, u1), p0.ov)
215	# N^4: O^3V^1 / N^4: O^3V^1
216	+ 2.0 * einsum("ia,jk->ijka",
217	einsum("il,la->ia",
218	einsum("lb,ib->il", u1, p0.ov), u1), d_oo)
219	# N^4: O^3V^1 / N^4: O^3V^1
220	+ 2.0 * einsum("ia,jk->ijka",
221	einsum("il,la->ia", einsum("ia,la->il", u1, u1),
222	p0.ov), d_oo)
223	).antisymmetrise(0, 1)
224	)
225	dm.oovv += (	2✔
226	# N^4: O^2V^2 / N^4: O^2V^2
227	+ 4.0 * einsum("ib,ja->ijab",
228	einsum("kc,ikbc->ib", u1, td2),
229	u1).antisymmetrise(0, 1).antisymmetrise(2, 3)
230	# N^5: O^3V^2 / N^4: O^2V^2
231	+ 2.0 * einsum("ijka,kb->ijab",
232	einsum("kc,ijac->ijka", u1, td2), u1).antisymmetrise(2, 3)
233	# N^5: O^3V^2 / N^4: O^2V^2
234	+ 2.0 * einsum("jk,ikab->ijab", einsum("jc,kc->jk", u1, u1), td2)
235	)
236	dm.ovov += (	2✔
237	# N^6: O^3V^3 / N^4: O^2V^2
238	- 4.0 * einsum("jkac,ikbc->iajb", u2, u2)
239	# N^5: O^2V^3 / N^4: O^2V^2
240	+ 2.0 * einsum("ab,ij->iajb", einsum("klac,klbc->ab", u2, u2), d_oo)
241	# N^6: O^4V^2 / N^4: O^2V^2
242	+ 1.0 * einsum("iklb,jkla->iajb", einsum("ld,ikbd->iklb", u1, t2),
243	einsum("lc,jkac->jkla", u1, t2))
244	# N^4: O^2V^2 / N^4: O^2V^2
245	+ 1.0 * einsum("ab,ij->iajb", einsum("ka,kb->ab", u1, u1), p0.oo)
246	# N^4: O^2V^2 / N^4: O^2V^2
247	- 1.0 * einsum("ij,ab->iajb", einsum("ic,jc->ij", u1, u1), p0.vv)
248	# N^6: O^3V^3 / N^4: O^2V^2
249	+ 1.0 * einsum("jlac,ilbc->iajb",
250	einsum("kl,jkac->jlac", einsum("kd,ld->kl", u1, u1), t2), t2)
251	# N^6: O^4V^2 / N^4: O^2V^2
252	+ 0.5 * einsum("ijla,lb->iajb",
253	einsum("ijkl,ka->ijla",
254	einsum("ikcd,jlcd->ijkl", t2, t2), u1), u1)
255	# N^6: O^4V^2 / N^4: O^2V^2
256	+ 0.5 * einsum("jklb,ikla->iajb",
257	einsum("jd,klbd->jklb", u1, t2),
258	einsum("ic,klac->ikla", u1, t2))
259	# N^4: O^2V^2 / N^4: O^2V^2
260	- 1.0 * einsum("ib,ja->iajb",
261	einsum("kd,ikbd->ib", u1, t2), einsum("lc,jlac->ja", u1, t2))
262	# N^5: O^2V^3 / N^4: O^2V^2
263	+ 1.0 * einsum("ab,ij->iajb",
264	einsum("lmad,lmbd->ab",
265	einsum("km,klad->lmad",
266	einsum("kc,mc->km", u1, u1), t2), t2), d_oo)
267	# N^4: O^2V^2 / N^4: O^2V^2
268	- 1.0 * einsum("ab,ij->iajb",
269	einsum("lb,la->ab",
270	einsum("kd,klbd->lb", u1, t2),
271	einsum("mc,lmac->la", u1, t2)), d_oo)
272	# N^5: O^3V^2 / N^4: O^2V^2
273	- 0.5 * einsum("ab,ij->iajb",
274	einsum("klmb,klma->ab",
275	einsum("md,klbd->klmb", u1, t2),
276	einsum("mc,klac->klma", u1, t2)), d_oo)
277	+ (
278	# N^4: O^2V^2 / N^4: O^2V^2
279	+ 1.0 * einsum("ib,ja->iajb", einsum("ic,bc->ib", u1, p0.vv), u1)
280	# N^4: O^2V^2 / N^4: O^2V^2
281	- 1.0 * einsum("ib,ja->iajb", einsum("kb,ik->ib", u1, p0.oo), u1)
282	# N^5: O^3V^2 / N^4: O^2V^2
283	+ 2.0 * einsum("ijlb,la->iajb",
284	einsum("jc,ilbc->ijlb",
285	einsum("kd,jkcd->jc", u1, t2), t2), u1)
286	# N^5: O^3V^2 / N^4: O^2V^2
287	+ 2.0 * einsum("kb,ijka->iajb",
288	einsum("ld,klbd->kb", u1, t2),
289	einsum("ic,jkac->ijka", u1, t2))
290	# N^6: O^4V^2 / N^4: O^2V^2
291	- 2.0 * einsum("ijkb,ka->iajb",
292	einsum("jklc,ilbc->ijkb",
293	einsum("kd,jlcd->jklc", u1, t2), t2), u1)
294	# N^6: O^4V^2 / N^4: O^2V^2
295	- 2.0 * einsum("ijlb,la->iajb",
296	einsum("ijkc,klbc->ijlb",
297	einsum("id,jkcd->ijkc", u1, t2), t2), u1)
298	# N^6: O^3V^3 / N^4: O^2V^2
299	- 2.0 * einsum("ikbc,jkac->iajb",
300	einsum("il,klbc->ikbc",
301	einsum("id,ld->il", u1, u1), t2), t2)
302	# N^4: O^2V^2 / N^4: O^2V^2
303	- 1.0 * einsum("ib,ja->iajb",
304	einsum("kc,ikbc->ib",
305	einsum("ld,klcd->kc", u1, t2), t2), u1)
306	# N^4: O^2V^2 / N^4: O^2V^2
307	- 1.0 * einsum("ab,ij->iajb",
308	einsum("kb,ka->ab",
309	einsum("kc,bc->kb", u1, p0.vv), u1), d_oo)
310	# N^4: O^2V^2 / N^4: O^2V^2
311	- 1.0 * einsum("ab,ij->iajb",
312	einsum("mb,ma->ab",
313	einsum("ld,lmbd->mb",
314	einsum("kc,klcd->ld", u1, t2),
315	t2), u1), d_oo)
316	).symmetrise([(0, 2), (1, 3)])
317	)
318	dm.ovvv += (	2✔
319	# N^6: O^3V^3 / N^4: O^1V^3
320	+ 1.0 * einsum("ijka,jkbc->iabc", einsum("id,jkad->ijka", u1, u2), t2)
321	# N^5: O^2V^3 / N^4: O^1V^3
322	+ 2.0 * einsum("ja,ijbc->iabc", einsum("kd,jkad->ja", u1, u2), t2)
323	+ (
324	# N^5: O^2V^3 / N^4: O^1V^3
325	+ 2.0 * einsum("ac,ib->iabc", einsum("jkad,jkcd->ac", u2, t2), u1)
326	# N^6: O^3V^3 / N^4: O^1V^3
327	- 4.0 * einsum("ikab,kc->iabc", einsum("jkad,ijbd->ikab", u2, t2), u1)
328	# N^4: O^1V^3 / N^4: O^1V^3
329	+ 2.0 * einsum("ac,ib->iabc", einsum("ja,jc->ac", u1, p0.ov), u1)
330	# N^4: O^1V^3 / N^4: O^1V^3
331	+ 2.0 * einsum("ac,ib->iabc", einsum("ja,jc->ac", u1, u1), p0.ov)
332	).antisymmetrise(2, 3)
333	)
334	dm.vvvv += (	2✔
335	# N^6: O^2V^4 / N^4: V^4
336	+ 2 * einsum("ijab,ijcd->abcd", u2, u2)
337	# N^6: O^2V^4 / N^4: V^4
338	- 1 * einsum("jkab,jkcd->abcd",
339	einsum("ij,ikab->jkab", einsum("ie,je->ij", u1, u1), t2), t2)
340	+ (
341	# N^4: V^4 / N^4: V^4
342	+ 4.0 * einsum("ac,bd->abcd", einsum("ia,ic->ac", u1, u1), p0.vv)
343	# N^6: O^3V^3 / N^4: V^4
344	+ 4.0 * einsum("jabc,jd->abcd",
345	einsum("ijbc,ia->jabc",
346	einsum("jkbe,ikce->ijbc", t2, t2), u1), u1)
347	).antisymmetrise(0, 1).antisymmetrise(2, 3)
348	+ (
349	# N^5: O^1V^4 / N^4: V^4
350	+ 4.0 * einsum("ka,kbcd->abcd",
351	einsum("ie,ikae->ka", u1, t2),
352	einsum("jb,jkcd->kbcd", u1, t2))
353	# N^6: O^3V^3 / N^4: V^4
354	+ 2.0 * einsum("ibcd,ia->abcd",
355	einsum("ijkb,jkcd->ibcd",
356	einsum("ie,jkbe->ijkb", u1, t2), t2), u1)
357	).antisymmetrise(0, 1).symmetrise([(0, 2), (1, 3)])
358	)
359	return dm	2✔
360
361
362	# dict controlling the dispatch of the state_diffdm function
363	DISPATCH = {	2✔
364	"isr0": diffdm_isr0_2p,
365	"isr1s": diffdm_isr1s_2p,
366	"isr1": diffdm_isr1_2p,
367	"isr2": diffdm_isr2_2p,
368	}
369
370
371	def state_diffdm_2p(method, ground_state, amplitude, intermediates=None):	2✔
372	"""
373	Compute the two-particle difference density matrix of an excited state
374	in the MO basis.
375
376	Parameters
377	----------
378	method : str, IsrMethod
379	The method to use for the computation (e.g. "isr2")
380	ground_state : LazyMp
381	The ground state upon which the excitation was based
382	amplitude : AmplitudeVector
383	The amplitude vector
384	intermediates : adcc.Intermediates
385	Intermediates from the ADC calculation to reuse
386	"""
387	if not isinstance(method, IsrMethod):	2!
NEW 388	method = IsrMethod(method)	×
389	if not isinstance(ground_state, LazyMp):	2!
390	raise TypeError("ground_state should be a LazyMp object.")	×
391	if not isinstance(amplitude, AmplitudeVector):	2!
392	raise TypeError("amplitude should be an AmplitudeVector object.")	×
393	if intermediates is None:	2✔
394	intermediates = Intermediates(ground_state)	2✔
395
396	if method.name not in DISPATCH:	2!
397	raise NotImplementedError("state_diffdm_2p is not implemented "	×
398	f"for {method.name}.")
399	else:
400	ret = DISPATCH[method.name](ground_state, amplitude, intermediates)	2✔
401	return ret.evaluate()	2✔

adc-connect / adcc / 24657644018

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous