adc-connect/adcc | Build 24441125763 | 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

    # ISR(1) 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

    # ISR(1) 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

        # ISR(1)-d

99

101

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

102

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

103

104

106

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

107

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

108

110

        # no doubles contribution

113

114

122

126

128

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

129

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

130

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

131

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

132

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

133

+ (

134

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

135

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

136

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

137

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

138

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

139

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

140

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

141

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

142

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

143

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

144

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

145

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

146

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

147

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

148

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

149

                           d_oo)

150

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

151

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

152

                           einsum("inbc,knbc->ik",

153

                                  einsum("mn,imbc->inbc",

154

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

155

                           d_oo)

156

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

157

+ (

158

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

159

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

160

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

161

                                  t2), u1)

162

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

163

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

164

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

165

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

166

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

167

+ (

168

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

169

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

170

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

171

                                                 p0.oo), d_oo)

172

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

173

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

174

                           einsum("kb,ib->ik",

175

                                  einsum("mc,kmbc->kb",

176

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

177

                           d_oo)

178

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

179

181

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

182

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

183

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

184

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

185

+ (

186

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

187

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

188

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

189

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

190

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

191

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

192

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

193

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

194

                           einsum("jlmb,lmab->ja",

195

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

196

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

197

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

198

                           einsum("lb,ilab->ia",

199

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

200

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

201

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

202

                           einsum("il,la->ia",

203

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

204

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

205

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

206

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

207

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

208

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

209

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

210

                           einsum("il,la->ia",

211

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

212

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

213

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

214

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

215

                                  p0.ov), d_oo)

216

        ).antisymmetrise(0, 1)

217

219

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

220

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

221

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

222

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

223

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

224

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

225

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

226

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

227

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

228

230

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

231

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

232

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

233

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

234

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

235

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

236

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

237

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

238

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

239

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

240

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

241

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

242

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

243

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

244

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

245

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

246

                       einsum("ijkl,ka->ijla",

247

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

248

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

249

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

250

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

251

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

252

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

253

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

254

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

255

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

256

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

257

                       einsum("lmad,lmbd->ab",

258

                              einsum("km,klad->lmad",

259

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

260

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

261

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

262

                       einsum("lb,la->ab",

263

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

264

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

265

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

266

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

267

                       einsum("klmb,klma->ab",

268

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

269

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

270

+ (

271

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

272

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

273

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

274

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

275

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

276

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

277

                           einsum("jc,ilbc->ijlb",

278

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

279

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

280

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

281

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

282

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

283

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

284

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

285

                           einsum("jklc,ilbc->ijkb",

286

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

287

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

288

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

289

                           einsum("ijkc,klbc->ijlb",

290

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

291

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

292

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

293

                           einsum("il,klbc->ikbc",

294

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

295

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

296

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

297

                           einsum("kc,ikbc->ib",

298

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

299

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

300

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

301

                           einsum("kb,ka->ab",

302

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

303

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

304

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

305

                           einsum("mb,ma->ab",

306

                                  einsum("ld,lmbd->mb",

307

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

308

                                         t2), u1), d_oo)

309

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

310

312

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

313

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

314

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

315

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

316

+ (

317

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

318

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

319

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

320

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

321

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

322

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

323

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

324

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

325

        ).antisymmetrise(2, 3)

326

328

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

329

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

330

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

331

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

332

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

333

+ (

334

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

335

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

336

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

337

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

338

                           einsum("ijbc,ia->jabc",

339

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

340

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

341

+ (

342

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

343

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

344

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

345

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

346

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

347

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

348

                           einsum("ijkb,jkcd->ibcd",

349

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

350

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

351

353

354

355

# dict controlling the dispatch of the state_diffdm function

357

    "isr0": diffdm_isr0_2p,

358

    "isr1": diffdm_isr1_2p,

359

    "isr2": diffdm_isr2_2p,

360

    "isr2x": diffdm_isr2_2p,      # same as ISR(2)

361

362

363

365

"""

366

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

367

    in the MO basis.

368

369

    Parameters

370

    ----------

371

    method : str, IsrMethod

372

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

373

    ground_state : LazyMp

374

        The ground state upon which the excitation was based

375

    amplitude : AmplitudeVector

376

        The amplitude vector

377

    intermediates : adcc.Intermediates

378

        Intermediates from the ADC calculation to reuse

379

"""

388

391

                                  f"for {method.name}.")

392

    else:

adc-connect / adcc / 24441125763

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	# ISR(1) 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_isr1_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	# ISR(1) 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	try:	2✔
96	# ISR(1)-d
97	check_doubles_amplitudes([b.o, b.o, b.v, b.v], amplitude)	2✔
98	u2 = amplitude.pphh	2✔
99
100	dm.ooov += (	2✔
101	# N^5: O^3V^2 / N^4: O^2V^2
102	- 2.0 * einsum("kb,ijab->ijka", u1, u2)
103	)
104
105	dm.ovvv += (	2✔
106	# N^5: O^2V^3 / N^4: O^1V^3
107	- 2.0 * einsum("ja,ijbc->iabc", u1, u2)
108	)
NEW 109	except ValueError:	×
110	# no doubles contribution
NEW 111	pass	×
112	return dm	2✔
113
114
115	def diffdm_isr2_2p(mp, amplitude, intermediates):	2✔
116	dm = diffdm_isr1_2p(mp, amplitude, intermediates) # Get ISR(1) result	2✔
117	check_doubles_amplitudes([b.o, b.o, b.v, b.v], amplitude)	2✔
118	u1, u2 = amplitude.ph, amplitude.pphh	2✔
119	hf = mp.reference_state	2✔
120	d_oo = zeros_like(hf.foo)	2✔
121	d_oo.set_mask("ii", 1)	2✔
122
123	t2 = mp.t2(b.oovv)	2✔
124	td2 = mp.td2(b.oovv)	2✔
125	p0 = mp.mp2_diffdm	2✔
126
127	dm.oooo += (	2✔
128	# N^6: O^4V^2 / N^4: O^2V^2
129	+ 2.0 * einsum("ijab,klab->ijkl", u2, u2)
130	# N^6: O^5V^1 / N^4: O^2V^2
131	- 1.0 * einsum("klmc,ijmc->ijkl", einsum("mb,klbc->klmc", u1, t2),
132	einsum("ma,ijac->ijmc", u1, t2))
133	+ (
134	# N^5: O^3V^2 / N^4: O^2V^2
135	- 8.0 * einsum("ik,jl->ijkl", einsum("imab,kmab->ik", u2, u2), d_oo)
136	# N^4: O^4 / N^4: O^4
137	+ 4.0 * einsum("il,jk->ijkl", einsum("ia,la->il", u1, u1), p0.oo)
138	# N^6: O^5V^1 / N^4: O^2V^2
139	+ 4.0 * einsum("jlmc,ikmc->ijkl", einsum("lb,jmbc->jlmc", u1, t2),
140	einsum("ia,kmac->ikmc", u1, t2))
141	# N^4: O^2V^2 / N^4: O^2V^2
142	+ 4.0 * einsum("il,jk->ijkl",
143	einsum("lc,ic->il", einsum("nb,lnbc->lc", u1, t2),
144	einsum("ma,imac->ic", u1, t2)), d_oo)
145	# N^5: O^3V^2 / N^4: O^2V^2
146	+ 4.0 * einsum("ik,jl->ijkl", einsum("kmnc,imnc->ik",
147	einsum("nb,kmbc->kmnc", u1, t2),
148	einsum("na,imac->imnc", u1, t2)),
149	d_oo)
150	# N^5: O^3V^2 / N^4: O^2V^2
151	+ 2.0 * einsum("ik,jl->ijkl",
152	einsum("inbc,knbc->ik",
153	einsum("mn,imbc->inbc",
154	einsum("ia,la->il", u1, u1), t2), t2),
155	d_oo)
156	).antisymmetrise(0, 1).antisymmetrise(2, 3)
157	+ (
158	# N^5: O^3V^2 / N^4: O^2V^2
159	+ 4.0 * einsum("iklb,jb->ijkl",
160	einsum("ic,klbc->iklb", einsum("ma,imac->ic", u1, t2),
161	t2), u1)
162	# N^6: O^4V^2 / N^4: O^2V^2
163	+ 2.0 * einsum("jklm,im->ijkl",
164	einsum("jmbc,klbc->jklm", t2, t2),
165	einsum("ia,la->il", u1, u1))
166	).antisymmetrise(0, 1).symmetrise([(0, 2), (1, 3)])
167	+ (
168	# N^4: O^4 / N^4: O^4
169	+ 4.0 * einsum("il,jk->ijkl", einsum("im,lm->il",
170	einsum("ia,la->il", u1, u1),
171	p0.oo), d_oo)
172	# N^4: O^2V^2 / N^4: O^2V^2
173	+ 4.0 * einsum("ik,jl->ijkl",
174	einsum("kb,ib->ik",
175	einsum("mc,kmbc->kb",
176	einsum("na,mnac->mc", u1, t2), t2), u1),
177	d_oo)
178	).antisymmetrise(0, 1).antisymmetrise(2, 3).symmetrise([(0, 2), (1, 3)])
179	)
180	dm.ooov += (	2✔
181	# N^6: O^4V^2 / N^4: O^2V^2
182	+ 1.0 * einsum("ijkl,la->ijka", einsum("klbc,ijbc->ijkl", u2, t2), u1)
183	# N^5: O^3V^2 / N^4: O^2V^2
184	+ 2.0 * einsum("kb,ijab->ijka", einsum("lc,klbc->kb", u1, u2), t2)
185	+ (
186	# N^4: O^2V^2 / N^4: O^2V^2
187	- 4.0 * einsum("ja,ik->ijka", einsum("lb,jlab->ja", u1, u2), d_oo)
188	# N^5: O^3V^2 / N^4: O^2V^2
189	+ 2.0 * einsum("jk,ia->ijka", einsum("klbc,jlbc->jk", u2, t2), u1)
190	# N^6: O^4V^2 / N^4: O^2V^2
191	- 4.0 * einsum("jklb,ilab->ijka", einsum("jc,klbc->jklb", u1, u2), t2)
192	# N^5: O^3V^2 / N^4: O^2V^2
193	+ 2.0 * einsum("ja,ik->ijka",
194	einsum("jlmb,lmab->ja",
195	einsum("jc,lmbc->jlmb", u1, u2), t2), d_oo)
196	# N^4: O^2V^2 / N^4: O^2V^2
197	- 4.0 * einsum("ia,jk->ijka",
198	einsum("lb,ilab->ia",
199	einsum("mc,lmbc->lb", u1, u2), t2), d_oo)
200	# N^5: O^3V^2 / N^4: O^2V^2
201	+ 2.0 * einsum("ia,jk->ijka",
202	einsum("il,la->ia",
203	einsum("lmbc,imbc->il", u2, t2), u1), d_oo)
204	# N^4: O^3V^1 / N^4: O^3V^1
205	+ 2.0 * einsum("jk,ia->ijka", einsum("kb,jb->jk", u1, p0.ov), u1)
206	# N^4: O^3V^1 / N^4: O^3V^1
207	+ 2.0 * einsum("jk,ia->ijka", einsum("jb,kb->jk", u1, u1), p0.ov)
208	# N^4: O^3V^1 / N^4: O^3V^1
209	+ 2.0 * einsum("ia,jk->ijka",
210	einsum("il,la->ia",
211	einsum("lb,ib->il", u1, p0.ov), u1), d_oo)
212	# N^4: O^3V^1 / N^4: O^3V^1
213	+ 2.0 * einsum("ia,jk->ijka",
214	einsum("il,la->ia", einsum("ia,la->il", u1, u1),
215	p0.ov), d_oo)
216	).antisymmetrise(0, 1)
217	)
218	dm.oovv += (	2✔
219	# N^4: O^2V^2 / N^4: O^2V^2
220	+ 4.0 * einsum("ib,ja->ijab",
221	einsum("kc,ikbc->ib", u1, td2),
222	u1).antisymmetrise(0, 1).antisymmetrise(2, 3)
223	# N^5: O^3V^2 / N^4: O^2V^2
224	+ 2.0 * einsum("ijka,kb->ijab",
225	einsum("kc,ijac->ijka", u1, td2), u1).antisymmetrise(2, 3)
226	# N^5: O^3V^2 / N^4: O^2V^2
227	+ 2.0 * einsum("jk,ikab->ijab", einsum("jc,kc->jk", u1, u1), td2)
228	)
229	dm.ovov += (	2✔
230	# N^6: O^3V^3 / N^4: O^2V^2
231	- 4.0 * einsum("jkac,ikbc->iajb", u2, u2)
232	# N^5: O^2V^3 / N^4: O^2V^2
233	+ 2.0 * einsum("ab,ij->iajb", einsum("klac,klbc->ab", u2, u2), d_oo)
234	# N^6: O^4V^2 / N^4: O^2V^2
235	+ 1.0 * einsum("iklb,jkla->iajb", einsum("ld,ikbd->iklb", u1, t2),
236	einsum("lc,jkac->jkla", u1, t2))
237	# N^4: O^2V^2 / N^4: O^2V^2
238	+ 1.0 * einsum("ab,ij->iajb", einsum("ka,kb->ab", u1, u1), p0.oo)
239	# N^4: O^2V^2 / N^4: O^2V^2
240	- 1.0 * einsum("ij,ab->iajb", einsum("ic,jc->ij", u1, u1), p0.vv)
241	# N^6: O^3V^3 / N^4: O^2V^2
242	+ 1.0 * einsum("jlac,ilbc->iajb",
243	einsum("kl,jkac->jlac", einsum("kd,ld->kl", u1, u1), t2), t2)
244	# N^6: O^4V^2 / N^4: O^2V^2
245	+ 0.5 * einsum("ijla,lb->iajb",
246	einsum("ijkl,ka->ijla",
247	einsum("ikcd,jlcd->ijkl", t2, t2), u1), u1)
248	# N^6: O^4V^2 / N^4: O^2V^2
249	+ 0.5 * einsum("jklb,ikla->iajb",
250	einsum("jd,klbd->jklb", u1, t2),
251	einsum("ic,klac->ikla", u1, t2))
252	# N^4: O^2V^2 / N^4: O^2V^2
253	- 1.0 * einsum("ib,ja->iajb",
254	einsum("kd,ikbd->ib", u1, t2), einsum("lc,jlac->ja", u1, t2))
255	# N^5: O^2V^3 / N^4: O^2V^2
256	+ 1.0 * einsum("ab,ij->iajb",
257	einsum("lmad,lmbd->ab",
258	einsum("km,klad->lmad",
259	einsum("kc,mc->km", u1, u1), t2), t2), d_oo)
260	# N^4: O^2V^2 / N^4: O^2V^2
261	- 1.0 * einsum("ab,ij->iajb",
262	einsum("lb,la->ab",
263	einsum("kd,klbd->lb", u1, t2),
264	einsum("mc,lmac->la", u1, t2)), d_oo)
265	# N^5: O^3V^2 / N^4: O^2V^2
266	- 0.5 * einsum("ab,ij->iajb",
267	einsum("klmb,klma->ab",
268	einsum("md,klbd->klmb", u1, t2),
269	einsum("mc,klac->klma", u1, t2)), d_oo)
270	+ (
271	# N^4: O^2V^2 / N^4: O^2V^2
272	+ 1.0 * einsum("ib,ja->iajb", einsum("ic,bc->ib", u1, p0.vv), u1)
273	# N^4: O^2V^2 / N^4: O^2V^2
274	- 1.0 * einsum("ib,ja->iajb", einsum("kb,ik->ib", u1, p0.oo), u1)
275	# N^5: O^3V^2 / N^4: O^2V^2
276	+ 2.0 * einsum("ijlb,la->iajb",
277	einsum("jc,ilbc->ijlb",
278	einsum("kd,jkcd->jc", u1, t2), t2), u1)
279	# N^5: O^3V^2 / N^4: O^2V^2
280	+ 2.0 * einsum("kb,ijka->iajb",
281	einsum("ld,klbd->kb", u1, t2),
282	einsum("ic,jkac->ijka", u1, t2))
283	# N^6: O^4V^2 / N^4: O^2V^2
284	- 2.0 * einsum("ijkb,ka->iajb",
285	einsum("jklc,ilbc->ijkb",
286	einsum("kd,jlcd->jklc", u1, t2), t2), u1)
287	# N^6: O^4V^2 / N^4: O^2V^2
288	- 2.0 * einsum("ijlb,la->iajb",
289	einsum("ijkc,klbc->ijlb",
290	einsum("id,jkcd->ijkc", u1, t2), t2), u1)
291	# N^6: O^3V^3 / N^4: O^2V^2
292	- 2.0 * einsum("ikbc,jkac->iajb",
293	einsum("il,klbc->ikbc",
294	einsum("id,ld->il", u1, u1), t2), t2)
295	# N^4: O^2V^2 / N^4: O^2V^2
296	- 1.0 * einsum("ib,ja->iajb",
297	einsum("kc,ikbc->ib",
298	einsum("ld,klcd->kc", u1, t2), t2), u1)
299	# N^4: O^2V^2 / N^4: O^2V^2
300	- 1.0 * einsum("ab,ij->iajb",
301	einsum("kb,ka->ab",
302	einsum("kc,bc->kb", u1, p0.vv), u1), d_oo)
303	# N^4: O^2V^2 / N^4: O^2V^2
304	- 1.0 * einsum("ab,ij->iajb",
305	einsum("mb,ma->ab",
306	einsum("ld,lmbd->mb",
307	einsum("kc,klcd->ld", u1, t2),
308	t2), u1), d_oo)
309	).symmetrise([(0, 2), (1, 3)])
310	)
311	dm.ovvv += (	2✔
312	# N^6: O^3V^3 / N^4: O^1V^3
313	+ 1.0 * einsum("ijka,jkbc->iabc", einsum("id,jkad->ijka", u1, u2), t2)
314	# N^5: O^2V^3 / N^4: O^1V^3
315	+ 2.0 * einsum("ja,ijbc->iabc", einsum("kd,jkad->ja", u1, u2), t2)
316	+ (
317	# N^5: O^2V^3 / N^4: O^1V^3
318	+ 2.0 * einsum("ac,ib->iabc", einsum("jkad,jkcd->ac", u2, t2), u1)
319	# N^6: O^3V^3 / N^4: O^1V^3
320	- 4.0 * einsum("ikab,kc->iabc", einsum("jkad,ijbd->ikab", u2, t2), u1)
321	# N^4: O^1V^3 / N^4: O^1V^3
322	+ 2.0 * einsum("ac,ib->iabc", einsum("ja,jc->ac", u1, p0.ov), u1)
323	# N^4: O^1V^3 / N^4: O^1V^3
324	+ 2.0 * einsum("ac,ib->iabc", einsum("ja,jc->ac", u1, u1), p0.ov)
325	).antisymmetrise(2, 3)
326	)
327	dm.vvvv += (	2✔
328	# N^6: O^2V^4 / N^4: V^4
329	+ 2 * einsum("ijab,ijcd->abcd", u2, u2)
330	# N^6: O^2V^4 / N^4: V^4
331	- 1 * einsum("jkab,jkcd->abcd",
332	einsum("ij,ikab->jkab", einsum("ie,je->ij", u1, u1), t2), t2)
333	+ (
334	# N^4: V^4 / N^4: V^4
335	+ 4.0 * einsum("ac,bd->abcd", einsum("ia,ic->ac", u1, u1), p0.vv)
336	# N^6: O^3V^3 / N^4: V^4
337	+ 4.0 * einsum("jabc,jd->abcd",
338	einsum("ijbc,ia->jabc",
339	einsum("jkbe,ikce->ijbc", t2, t2), u1), u1)
340	).antisymmetrise(0, 1).antisymmetrise(2, 3)
341	+ (
342	# N^5: O^1V^4 / N^4: V^4
343	+ 4.0 * einsum("ka,kbcd->abcd",
344	einsum("ie,ikae->ka", u1, t2),
345	einsum("jb,jkcd->kbcd", u1, t2))
346	# N^6: O^3V^3 / N^4: V^4
347	+ 2.0 * einsum("ibcd,ia->abcd",
348	einsum("ijkb,jkcd->ibcd",
349	einsum("ie,jkbe->ijkb", u1, t2), t2), u1)
350	).antisymmetrise(0, 1).symmetrise([(0, 2), (1, 3)])
351	)
352	return dm	2✔
353
354
355	# dict controlling the dispatch of the state_diffdm function
356	DISPATCH = {	2✔
357	"isr0": diffdm_isr0_2p,
358	"isr1": diffdm_isr1_2p,
359	"isr2": diffdm_isr2_2p,
360	"isr2x": diffdm_isr2_2p, # same as ISR(2)
361	}
362
363
364	def state_diffdm_2p(method, ground_state, amplitude, intermediates=None):	2✔
365	"""
366	Compute the two-particle difference density matrix of an excited state
367	in the MO basis.
368
369	Parameters
370	----------
371	method : str, IsrMethod
372	The method to use for the computation (e.g. "isr2")
373	ground_state : LazyMp
374	The ground state upon which the excitation was based
375	amplitude : AmplitudeVector
376	The amplitude vector
377	intermediates : adcc.Intermediates
378	Intermediates from the ADC calculation to reuse
379	"""
380	if not isinstance(method, IsrMethod):	2!
NEW UNCOV 381	method = IsrMethod(method)	×
382	if not isinstance(ground_state, LazyMp):	2!
UNCOV 383	raise TypeError("ground_state should be a LazyMp object.")	×
384	if not isinstance(amplitude, AmplitudeVector):	2!
UNCOV 385	raise TypeError("amplitude should be an AmplitudeVector object.")	×
386	if intermediates is None:	2✔
387	intermediates = Intermediates(ground_state)	2✔
388
389	if method.name not in DISPATCH:	2!
390	raise NotImplementedError("state_diffdm_2p is not implemented "	×
391	f"for {method.name}.")
392	else:
393	ret = DISPATCH[method.name](ground_state, amplitude, intermediates)	2✔
394	return ret.evaluate()	2✔

adc-connect / adcc / 24441125763

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