| 1 | //
|
|---|
| 2 | // compute_vxb_a_asymm.cc
|
|---|
| 3 | //
|
|---|
| 4 | // Copyright (C) 2004 Edward Valeev
|
|---|
| 5 | //
|
|---|
| 6 | // Author: Edward Valeev <edward.valeev@chemistry.gatech.edu>
|
|---|
| 7 | // Maintainer: EV
|
|---|
| 8 | //
|
|---|
| 9 | // This file is part of the SC Toolkit.
|
|---|
| 10 | //
|
|---|
| 11 | // The SC Toolkit is free software; you can redistribute it and/or modify
|
|---|
| 12 | // it under the terms of the GNU Library General Public License as published by
|
|---|
| 13 | // the Free Software Foundation; either version 2, or (at your option)
|
|---|
| 14 | // any later version.
|
|---|
| 15 | //
|
|---|
| 16 | // The SC Toolkit is distributed in the hope that it will be useful,
|
|---|
| 17 | // but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|---|
| 18 | // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|---|
| 19 | // GNU Library General Public License for more details.
|
|---|
| 20 | //
|
|---|
| 21 | // You should have received a copy of the GNU Library General Public License
|
|---|
| 22 | // along with the SC Toolkit; see the file COPYING.LIB. If not, write to
|
|---|
| 23 | // the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
|
|---|
| 24 | //
|
|---|
| 25 | // The U.S. Government is granted a limited license as per AL 91-7.
|
|---|
| 26 | //
|
|---|
| 27 |
|
|---|
| 28 | #include <stdexcept>
|
|---|
| 29 | #include <sstream>
|
|---|
| 30 | #include <stdlib.h>
|
|---|
| 31 | #include <math.h>
|
|---|
| 32 | #include <limits.h>
|
|---|
| 33 |
|
|---|
| 34 | #include <scconfig.h>
|
|---|
| 35 | #include <util/misc/formio.h>
|
|---|
| 36 | #include <util/misc/timer.h>
|
|---|
| 37 | #include <util/class/class.h>
|
|---|
| 38 | #include <util/state/state.h>
|
|---|
| 39 | #include <util/state/state_text.h>
|
|---|
| 40 | #include <util/state/state_bin.h>
|
|---|
| 41 | #include <math/scmat/matrix.h>
|
|---|
| 42 | #include <chemistry/molecule/molecule.h>
|
|---|
| 43 | #include <chemistry/qc/basis/integral.h>
|
|---|
| 44 | #include <chemistry/qc/mbpt/bzerofast.h>
|
|---|
| 45 | #include <chemistry/qc/mbptr12/r12ia.h>
|
|---|
| 46 | #include <chemistry/qc/mbptr12/vxb_eval_info.h>
|
|---|
| 47 | #include <chemistry/qc/mbptr12/pairiter.h>
|
|---|
| 48 | #include <chemistry/qc/mbptr12/r12int_eval.h>
|
|---|
| 49 |
|
|---|
| 50 | using namespace std;
|
|---|
| 51 | using namespace sc;
|
|---|
| 52 |
|
|---|
| 53 | #define PRINT_R12_INTERMED 0
|
|---|
| 54 |
|
|---|
| 55 | void
|
|---|
| 56 | R12IntEval::contrib_to_VXB_a_asymm_(const std::string& tform_name)
|
|---|
| 57 | {
|
|---|
| 58 | if (evaluated_)
|
|---|
| 59 | return;
|
|---|
| 60 | LinearR12::ABSMethod abs_method = r12info_->abs_method();
|
|---|
| 61 | Ref<MessageGrp> msg = r12info_->msg();
|
|---|
| 62 | Ref<MemoryGrp> mem = r12info_->mem();
|
|---|
| 63 | Ref<ThreadGrp> thr = r12info_->thr();
|
|---|
| 64 | const int num_te_types = 4;
|
|---|
| 65 | enum te_types {eri=0, r12=1, r12t1=2, r12t2=3};
|
|---|
| 66 |
|
|---|
| 67 | tim_enter("mp2-r12a intermeds (asymmetric term)");
|
|---|
| 68 |
|
|---|
| 69 | int me = msg->me();
|
|---|
| 70 | int nproc = msg->n();
|
|---|
| 71 |
|
|---|
| 72 | // Do the AO->MO transform
|
|---|
| 73 | Ref<TwoBodyMOIntsTransform> ikjy_tform = get_tform_(tform_name);
|
|---|
| 74 | Ref<R12IntsAcc> ijky_acc = ikjy_tform->ints_acc();
|
|---|
| 75 | if (ijky_acc.null() || !ijky_acc->is_committed())
|
|---|
| 76 | ikjy_tform->compute();
|
|---|
| 77 | if (!ijky_acc->is_active())
|
|---|
| 78 | ijky_acc->activate();
|
|---|
| 79 | if (num_te_types != ijky_acc->num_te_types())
|
|---|
| 80 | throw std::runtime_error("R12IntEval::contrib_to_VXB_a_asymm_() -- number of MO integral types is wrong");
|
|---|
| 81 |
|
|---|
| 82 | Ref<MOIndexSpace> mospace1 = ikjy_tform->space2();
|
|---|
| 83 | Ref<MOIndexSpace> mospace2 = ikjy_tform->space4();
|
|---|
| 84 |
|
|---|
| 85 | ostringstream oss;
|
|---|
| 86 | oss << "\"" << mospace1->name() << "\"/\"" << mospace2->name() << "\"";
|
|---|
| 87 | std::string label = oss.str();
|
|---|
| 88 | ExEnv::out0() << endl << indent
|
|---|
| 89 | << "Entered " << label
|
|---|
| 90 | << " A (GEBC) intermediates evaluator" << endl;
|
|---|
| 91 | ExEnv::out0() << incindent;
|
|---|
| 92 |
|
|---|
| 93 | const int rank2 = mospace1->rank();
|
|---|
| 94 | const int rank4 = mospace2->rank();
|
|---|
| 95 |
|
|---|
| 96 | /*--------------------------------
|
|---|
| 97 | Compute MP2-R12/A intermediates
|
|---|
| 98 | and collect on node0
|
|---|
| 99 | --------------------------------*/
|
|---|
| 100 | ExEnv::out0() << indent << "Begin computation of intermediates" << endl;
|
|---|
| 101 | tim_enter("intermediates");
|
|---|
| 102 | SpatialMOPairIter_eq ij_iter(r12info_->act_occ_space());
|
|---|
| 103 | SpatialMOPairIter_eq kl_iter(r12info_->act_occ_space());
|
|---|
| 104 | int naa = ij_iter.nij_aa(); // Number of alpha-alpha pairs (i > j)
|
|---|
| 105 | int nab = ij_iter.nij_ab(); // Number of alpha-beta pairs
|
|---|
| 106 | if (debug_) {
|
|---|
| 107 | ExEnv::out0() << indent << "naa = " << naa << endl;
|
|---|
| 108 | ExEnv::out0() << indent << "nab = " << nab << endl;
|
|---|
| 109 | }
|
|---|
| 110 |
|
|---|
| 111 | // Compute intermediates
|
|---|
| 112 | if (debug_)
|
|---|
| 113 | ExEnv::out0() << indent << "Ready to compute MP2-R12/A (GEBC) intermediates" << endl;
|
|---|
| 114 |
|
|---|
| 115 | // Compute the number of tasks that have full access to the integrals
|
|---|
| 116 | // and split the work among them
|
|---|
| 117 | vector<int> proc_with_ints;
|
|---|
| 118 | int nproc_with_ints = tasks_with_ints_(ijky_acc,proc_with_ints);
|
|---|
| 119 |
|
|---|
| 120 |
|
|---|
| 121 | //////////////////////////////////////////////////////////////
|
|---|
| 122 | //
|
|---|
| 123 | // Evaluation of the intermediates proceeds as follows:
|
|---|
| 124 | //
|
|---|
| 125 | // loop over batches of kl, k >= l, 0<=k,l<nocc_act
|
|---|
| 126 | // load (kl|xy), (kl| [T1,r12] |xy), and (lk| [T1,r12] |xy)
|
|---|
| 127 | // (aka kl-sets) into memory
|
|---|
| 128 | //
|
|---|
| 129 | // loop over batches of ij, i>=j, 0<=i,j<nocc_act
|
|---|
| 130 | // load (ij|r12|xy) into memory
|
|---|
| 131 | // (aka ij-sets) into memory
|
|---|
| 132 | // compute V[ij][kl] and T[ij][kl] for all ij and kl in
|
|---|
| 133 | // the "direct product" batch
|
|---|
| 134 | // end ij loop
|
|---|
| 135 | // end kl loop
|
|---|
| 136 | //
|
|---|
| 137 | /////////////////////////////////////////////////////////////////////////////////
|
|---|
| 138 |
|
|---|
| 139 | if (ijky_acc->has_access(me)) {
|
|---|
| 140 |
|
|---|
| 141 | for(kl_iter.start();int(kl_iter);kl_iter.next()) {
|
|---|
| 142 |
|
|---|
| 143 | const int kl = kl_iter.ij();
|
|---|
| 144 | // Figure out if this task will handle this kl
|
|---|
| 145 | int kl_proc = kl%nproc_with_ints;
|
|---|
| 146 | if (kl_proc != proc_with_ints[me])
|
|---|
| 147 | continue;
|
|---|
| 148 | const int k = kl_iter.i();
|
|---|
| 149 | const int l = kl_iter.j();
|
|---|
| 150 | const int kl_aa = kl_iter.ij_aa();
|
|---|
| 151 | const int kl_ab = kl_iter.ij_ab();
|
|---|
| 152 | const int lk_ab = kl_iter.ij_ba();
|
|---|
| 153 |
|
|---|
| 154 | if (debug_)
|
|---|
| 155 | ExEnv::outn() << indent << "task " << me << ": working on (k,l) = " << k << "," << l << " " << endl;
|
|---|
| 156 |
|
|---|
| 157 | // Get (|1/r12|), (|r12|), (|[r12,T1]|), and (|[r12,T2]|) integrals
|
|---|
| 158 | tim_enter("MO ints retrieve");
|
|---|
| 159 | double *klox_buf_eri = ijky_acc->retrieve_pair_block(k,l,R12IntsAcc::eri);
|
|---|
| 160 | double *klox_buf_r12 = ijky_acc->retrieve_pair_block(k,l,R12IntsAcc::r12);
|
|---|
| 161 | double *klox_buf_r12t1 = ijky_acc->retrieve_pair_block(k,l,R12IntsAcc::r12t1);
|
|---|
| 162 | double *klox_buf_r12t2 = ijky_acc->retrieve_pair_block(k,l,R12IntsAcc::r12t2);
|
|---|
| 163 |
|
|---|
| 164 | double *lkox_buf_eri = ijky_acc->retrieve_pair_block(l,k,R12IntsAcc::eri);
|
|---|
| 165 | double *lkox_buf_r12 = ijky_acc->retrieve_pair_block(l,k,R12IntsAcc::r12);
|
|---|
| 166 | double *lkox_buf_r12t1 = ijky_acc->retrieve_pair_block(l,k,R12IntsAcc::r12t1);
|
|---|
| 167 | double *lkox_buf_r12t2 = ijky_acc->retrieve_pair_block(l,k,R12IntsAcc::r12t2);
|
|---|
| 168 | tim_exit("MO ints retrieve");
|
|---|
| 169 |
|
|---|
| 170 | if (debug_)
|
|---|
| 171 | ExEnv::outn() << indent << "task " << me << ": obtained kl blocks" << endl;
|
|---|
| 172 |
|
|---|
| 173 | // to avoid every task hitting same ij at the same time, stagger ij-accesses, i.e. each kl task will start with ij=kl+1
|
|---|
| 174 | for(ij_iter.start(kl+1);int(ij_iter);ij_iter.next()) {
|
|---|
| 175 |
|
|---|
| 176 | const int i = ij_iter.i();
|
|---|
| 177 | const int j = ij_iter.j();
|
|---|
| 178 | const int ij_aa = ij_iter.ij_aa();
|
|---|
| 179 | const int ij_ab = ij_iter.ij_ab();
|
|---|
| 180 | const int ji_ab = ij_iter.ij_ba();
|
|---|
| 181 |
|
|---|
| 182 | if (debug_)
|
|---|
| 183 | ExEnv::outn() << indent << "task " << me << ": (k,l) = " << k << "," << l << ": (i,j) = " << i << "," << j << endl;
|
|---|
| 184 |
|
|---|
| 185 | tim_enter("MO ints retrieve");
|
|---|
| 186 | double *ijox_buf_r12 = ijky_acc->retrieve_pair_block(i,j,R12IntsAcc::r12);
|
|---|
| 187 | double *jiox_buf_r12 = ijky_acc->retrieve_pair_block(j,i,R12IntsAcc::r12);
|
|---|
| 188 | tim_exit("MO ints retrieve");
|
|---|
| 189 |
|
|---|
| 190 | if (debug_)
|
|---|
| 191 | ExEnv::outn() << indent << "task " << me << ": obtained ij blocks" << endl;
|
|---|
| 192 |
|
|---|
| 193 |
|
|---|
| 194 | tim_enter("MO ints contraction");
|
|---|
| 195 | double Vaa_ijkl, Vab_ijkl, Vab_jikl, Vab_ijlk, Vab_jilk;
|
|---|
| 196 | double Xaa_ijkl, Xab_ijkl, Xab_jikl, Xab_ijlk, Xab_jilk;
|
|---|
| 197 | double Taa_ijkl, Tab_ijkl, Tab_jikl, Tab_ijlk, Tab_jilk;
|
|---|
| 198 | Vaa_ijkl = Vab_ijkl = Vab_jikl = Vab_ijlk = Vab_jilk = 0.0;
|
|---|
| 199 | Xaa_ijkl = Xab_ijkl = Xab_jikl = Xab_ijlk = Xab_jilk = 0.0;
|
|---|
| 200 | Taa_ijkl = Tab_ijkl = Tab_jikl = Tab_ijlk = Tab_jilk = 0.0;
|
|---|
| 201 | for(int o=0; o<rank2; o++) {
|
|---|
| 202 | const double pfac_xy = 1.0;
|
|---|
| 203 | for(int x=0; x<rank4; x++) {
|
|---|
| 204 | int ox_offset = o*rank4 + x;
|
|---|
| 205 | double ij_r12_ox = ijox_buf_r12[ox_offset];
|
|---|
| 206 | double ji_r12_ox = jiox_buf_r12[ox_offset];
|
|---|
| 207 | double kl_eri_ox = klox_buf_eri[ox_offset];
|
|---|
| 208 | double lk_eri_ox = lkox_buf_eri[ox_offset];
|
|---|
| 209 | Vab_ijkl -= pfac_xy * (ij_r12_ox * kl_eri_ox + ji_r12_ox * lk_eri_ox);
|
|---|
| 210 | if (ij_ab != ji_ab)
|
|---|
| 211 | Vab_jikl -= pfac_xy * (ji_r12_ox * kl_eri_ox + ij_r12_ox * lk_eri_ox);
|
|---|
| 212 | if (kl_ab != lk_ab)
|
|---|
| 213 | Vab_ijlk -= pfac_xy * (ij_r12_ox * lk_eri_ox + ji_r12_ox * kl_eri_ox);
|
|---|
| 214 | if (ij_ab != ji_ab && kl_ab != lk_ab) {
|
|---|
| 215 | Vab_jilk -= pfac_xy * (ij_r12_ox * kl_eri_ox + ji_r12_ox * lk_eri_ox);
|
|---|
| 216 | }
|
|---|
| 217 | if (ij_aa != -1 && kl_aa != -1) {
|
|---|
| 218 | Vaa_ijkl -= pfac_xy * (ij_r12_ox - ji_r12_ox)*(kl_eri_ox - lk_eri_ox);
|
|---|
| 219 | }
|
|---|
| 220 | double kl_r12_ox = klox_buf_r12[ox_offset];
|
|---|
| 221 | double lk_r12_ox = lkox_buf_r12[ox_offset];
|
|---|
| 222 | Xab_ijkl -= pfac_xy * (ij_r12_ox * kl_r12_ox + ji_r12_ox * lk_r12_ox);
|
|---|
| 223 | if (ij_ab != ji_ab)
|
|---|
| 224 | Xab_jikl -= pfac_xy * (ji_r12_ox * kl_r12_ox + ij_r12_ox * lk_r12_ox);
|
|---|
| 225 | if (kl_ab != lk_ab)
|
|---|
| 226 | Xab_ijlk -= pfac_xy * (ij_r12_ox * lk_r12_ox + ji_r12_ox * kl_r12_ox);
|
|---|
| 227 | if (ij_ab != ji_ab && kl_ab != lk_ab) {
|
|---|
| 228 | Xab_jilk -= pfac_xy * (ij_r12_ox * kl_r12_ox + ji_r12_ox * lk_r12_ox);
|
|---|
| 229 | }
|
|---|
| 230 | if (ij_aa != -1 && kl_aa != -1) {
|
|---|
| 231 | Xaa_ijkl -= pfac_xy * (ij_r12_ox - ji_r12_ox)*(kl_r12_ox - lk_r12_ox);
|
|---|
| 232 | }
|
|---|
| 233 | double kl_r12t1_ox = klox_buf_r12t1[ox_offset];
|
|---|
| 234 | double kl_r12t2_ox = klox_buf_r12t2[ox_offset];
|
|---|
| 235 | double lk_r12t1_ox = lkox_buf_r12t1[ox_offset];
|
|---|
| 236 | double lk_r12t2_ox = lkox_buf_r12t2[ox_offset];
|
|---|
| 237 | double kl_Tr12_ox = -kl_r12t1_ox-kl_r12t2_ox;
|
|---|
| 238 | double lk_Tr12_ox = -lk_r12t1_ox-lk_r12t2_ox;
|
|---|
| 239 | Tab_ijkl += pfac_xy * (ij_r12_ox * kl_Tr12_ox + ji_r12_ox * lk_Tr12_ox);
|
|---|
| 240 | if (ij_ab != ji_ab)
|
|---|
| 241 | Tab_jikl += pfac_xy * (ji_r12_ox * kl_Tr12_ox + ij_r12_ox * lk_Tr12_ox);
|
|---|
| 242 | if (kl_ab != lk_ab)
|
|---|
| 243 | Tab_ijlk += pfac_xy * (ij_r12_ox * lk_Tr12_ox + ji_r12_ox * kl_Tr12_ox);
|
|---|
| 244 | if (ij_ab != ji_ab && kl_ab != lk_ab) {
|
|---|
| 245 | Tab_jilk += pfac_xy * (ij_r12_ox * kl_Tr12_ox + ji_r12_ox * lk_Tr12_ox);
|
|---|
| 246 | }
|
|---|
| 247 | if (ij_aa != -1 && kl_aa != -1) {
|
|---|
| 248 | Taa_ijkl += pfac_xy * (ij_r12_ox - ji_r12_ox)*(kl_Tr12_ox - lk_Tr12_ox);
|
|---|
| 249 | }
|
|---|
| 250 | }
|
|---|
| 251 | }
|
|---|
| 252 | Vab_.accumulate_element(ij_ab,kl_ab,Vab_ijkl);
|
|---|
| 253 | if (ij_ab != ji_ab)
|
|---|
| 254 | Vab_.accumulate_element(ji_ab,kl_ab,Vab_jikl);
|
|---|
| 255 | if (kl_ab != lk_ab)
|
|---|
| 256 | Vab_.accumulate_element(ij_ab,lk_ab,Vab_ijlk);
|
|---|
| 257 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 258 | Vab_.accumulate_element(ji_ab,lk_ab,Vab_jilk);
|
|---|
| 259 | if (ij_aa != -1 && kl_aa != -1)
|
|---|
| 260 | Vaa_.accumulate_element(ij_aa,kl_aa,Vaa_ijkl);
|
|---|
| 261 | Xab_.accumulate_element(ij_ab,kl_ab,Xab_ijkl);
|
|---|
| 262 | if (ij_ab != ji_ab)
|
|---|
| 263 | Xab_.accumulate_element(ji_ab,kl_ab,Xab_jikl);
|
|---|
| 264 | if (kl_ab != lk_ab)
|
|---|
| 265 | Xab_.accumulate_element(ij_ab,lk_ab,Xab_ijlk);
|
|---|
| 266 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 267 | Xab_.accumulate_element(ji_ab,lk_ab,Xab_jilk);
|
|---|
| 268 | if (ij_aa != -1 && kl_aa != -1)
|
|---|
| 269 | Xaa_.accumulate_element(ij_aa,kl_aa,Xaa_ijkl);
|
|---|
| 270 | Bab_.accumulate_element(ij_ab,kl_ab,Tab_ijkl);
|
|---|
| 271 | if (ij_ab != ji_ab)
|
|---|
| 272 | Bab_.accumulate_element(ji_ab,kl_ab,Tab_jikl);
|
|---|
| 273 | if (kl_ab != lk_ab)
|
|---|
| 274 | Bab_.accumulate_element(ij_ab,lk_ab,Tab_ijlk);
|
|---|
| 275 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 276 | Bab_.accumulate_element(ji_ab,lk_ab,Tab_jilk);
|
|---|
| 277 | if (ij_aa != -1 && kl_aa != -1)
|
|---|
| 278 | Baa_.accumulate_element(ij_aa,kl_aa,Taa_ijkl);
|
|---|
| 279 | tim_exit("MO ints contraction");
|
|---|
| 280 |
|
|---|
| 281 | #if PRINT_R12_INTERMED
|
|---|
| 282 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 283 | printf("Vaa[%d][%d] = %lf\n",ij_aa,kl_aa,Vaa_ij[kl_aa]);
|
|---|
| 284 | printf("Vab[%d][%d] = %lf\n",ij_ab,kl_ab,Vab_ij[kl_ab]);
|
|---|
| 285 | if (ij_ab != ji_ab)
|
|---|
| 286 | printf("Vab[%d][%d] = %lf\n",ji_ab,kl_ab,Vab_ji[kl_ab]);
|
|---|
| 287 | if (kl_ab != lk_ab)
|
|---|
| 288 | printf("Vab[%d][%d] = %lf\n",ij_ab,lk_ab,Vab_ij[lk_ab]);
|
|---|
| 289 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 290 | printf("Vab[%d][%d] = %lf\n",ji_ab,lk_ab,Vab_ji[lk_ab]);
|
|---|
| 291 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 292 | printf("Xaa[%d][%d] = %lf\n",ij_aa,kl_aa,Xaa_ij[kl_aa]);
|
|---|
| 293 | printf("Xab[%d][%d] = %lf\n",ij_ab,kl_ab,Xab_ij[kl_ab]);
|
|---|
| 294 | if (ij_ab != ji_ab)
|
|---|
| 295 | printf("Xab[%d][%d] = %lf\n",ji_ab,kl_ab,Xab_ji[kl_ab]);
|
|---|
| 296 | if (kl_ab != lk_ab)
|
|---|
| 297 | printf("Xab[%d][%d] = %lf\n",ij_ab,lk_ab,Xab_ij[lk_ab]);
|
|---|
| 298 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 299 | printf("Xab[%d][%d] = %lf\n",ji_ab,lk_ab,Xab_ji[lk_ab]);
|
|---|
| 300 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 301 | printf("Taa[%d][%d] = %lf\n",ij_aa,kl_aa,Taa_ij[kl_aa]);
|
|---|
| 302 | printf("Tab[%d][%d] = %lf\n",ij_ab,kl_ab,Tab_ij[kl_ab]);
|
|---|
| 303 | if (ij_ab != ji_ab)
|
|---|
| 304 | printf("Tab[%d][%d] = %lf\n",ji_ab,kl_ab,Tab_ji[kl_ab]);
|
|---|
| 305 | if (kl_ab != lk_ab)
|
|---|
| 306 | printf("Tab[%d][%d] = %lf\n",ij_ab,lk_ab,Tab_ij[lk_ab]);
|
|---|
| 307 | if (ij_ab != ji_ab && kl_ab != lk_ab)
|
|---|
| 308 | printf("Tab[%d][%d] = %lf\n",ji_ab,lk_ab,Tab_ji[lk_ab]);
|
|---|
| 309 | #endif
|
|---|
| 310 | ijky_acc->release_pair_block(i,j,R12IntsAcc::r12);
|
|---|
| 311 | ijky_acc->release_pair_block(j,i,R12IntsAcc::r12);
|
|---|
| 312 | }
|
|---|
| 313 | ijky_acc->release_pair_block(k,l,R12IntsAcc::eri);
|
|---|
| 314 | ijky_acc->release_pair_block(k,l,R12IntsAcc::r12);
|
|---|
| 315 | ijky_acc->release_pair_block(k,l,R12IntsAcc::r12t1);
|
|---|
| 316 | ijky_acc->release_pair_block(k,l,R12IntsAcc::r12t2);
|
|---|
| 317 | ijky_acc->release_pair_block(l,k,R12IntsAcc::eri);
|
|---|
| 318 | ijky_acc->release_pair_block(l,k,R12IntsAcc::r12);
|
|---|
| 319 | ijky_acc->release_pair_block(l,k,R12IntsAcc::r12t1);
|
|---|
| 320 | ijky_acc->release_pair_block(l,k,R12IntsAcc::r12t2);
|
|---|
| 321 | }
|
|---|
| 322 | }
|
|---|
| 323 | // Tasks that don't do any work here still need to create these timers
|
|---|
| 324 | tim_enter("MO ints retrieve");
|
|---|
| 325 | tim_exit("MO ints retrieve");
|
|---|
| 326 | tim_enter("MO ints contraction");
|
|---|
| 327 | tim_exit("MO ints contraction");
|
|---|
| 328 |
|
|---|
| 329 | tim_exit("intermediates");
|
|---|
| 330 | ExEnv::out0() << indent << "End of computation of intermediates" << endl;
|
|---|
| 331 | ijky_acc->deactivate();
|
|---|
| 332 |
|
|---|
| 333 | // Symmetrize B intermediate
|
|---|
| 334 | for(int ij=0;ij<naa;ij++)
|
|---|
| 335 | for(int kl=0;kl<=ij;kl++) {
|
|---|
| 336 | double belem = 0.5*(Baa_->get_element(ij,kl) + Baa_->get_element(kl,ij));
|
|---|
| 337 | Baa_->set_element(ij,kl,belem);
|
|---|
| 338 | Baa_->set_element(kl,ij,belem);
|
|---|
| 339 | }
|
|---|
| 340 |
|
|---|
| 341 | for(int ij=0;ij<nab;ij++)
|
|---|
| 342 | for(int kl=0;kl<=ij;kl++) {
|
|---|
| 343 | double belem = 0.5*(Bab_->get_element(ij,kl) + Bab_->get_element(kl,ij));
|
|---|
| 344 | Bab_->set_element(ij,kl,belem);
|
|---|
| 345 | Bab_->set_element(kl,ij,belem);
|
|---|
| 346 | }
|
|---|
| 347 |
|
|---|
| 348 | globally_sum_intermeds_();
|
|---|
| 349 |
|
|---|
| 350 | ExEnv::out0() << decindent;
|
|---|
| 351 | ExEnv::out0() << indent
|
|---|
| 352 | << "Exited " << label
|
|---|
| 353 | << " A (GEBC) intermediates evaluator" << endl;
|
|---|
| 354 |
|
|---|
| 355 | tim_exit("mp2-r12a intermeds (asymmetric term)");
|
|---|
| 356 | checkpoint_();
|
|---|
| 357 | }
|
|---|
| 358 |
|
|---|
| 359 | ////////////////////////////////////////////////////////////////////////////
|
|---|
| 360 |
|
|---|
| 361 | // Local Variables:
|
|---|
| 362 | // mode: c++
|
|---|
| 363 | // c-file-style: "CLJ-CONDENSED"
|
|---|
| 364 | // End:
|
|---|
