- Main Page
- Related Pages
- Classes
- Files

# lbgbuild.h

    00001 //
    00002 // lbgbuild.h --- definitino of the load-balanced local G matrix builder
    00003 //
    00004 // Copyright (C) 1996 Limit Point Systems, Inc.
    00005 //
    00006 // Author: Edward Seidl <seidl@janed.com>
    00007 // Maintainer: LPS
    00008 //
    00009 // This file is part of the SC Toolkit.
    00010 //
    00011 // The SC Toolkit is free software; you can redistribute it and/or modify
    00012 // it under the terms of the GNU Library General Public License as published by
    00013 // the Free Software Foundation; either version 2, or (at your option)
    00014 // any later version.
    00015 //
    00016 // The SC Toolkit is distributed in the hope that it will be useful,
    00017 // but WITHOUT ANY WARRANTY; without even the implied warranty of
    00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    00019 // GNU Library General Public License for more details.
    00020 //
    00021 // You should have received a copy of the GNU Library General Public License
    00022 // along with the SC Toolkit; see the file COPYING.LIB.  If not, write to
    00023 // the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
    00024 //
    00025 // The U.S. Government is granted a limited license as per AL 91-7.
    00026 //
    00027
    00028 #ifndef _chemistry_qc_scf_lbgbuild_h
    00029 #define _chemistry_qc_scf_lbgbuild_h
    00030
    00031 #ifdef __GNUC__
    00032 #pragma interface
    00033 #endif
    00034
    00035 #include <chemistry/qc/scf/gbuild.h>
    00036
    00037 namespace sc {
    00038
    00039 template<class T>
    00040 class LocalLBGBuild : public GBuild<T> {
    00041   protected:
    00042     Ref<MessageGrp> grp_;
    00043     Ref<TwoBodyInt> tbi_;
    00044     Ref<Integral> integral_;
    00045     Ref<GaussianBasisSet> gbs_;
    00046     signed char *pmax;
    00047
    00048   public:
    00049     LocalLBGBuild(T& t, const Ref<TwoBodyInt>& tbi, const Ref<Integral>& ints,
    00050                 const Ref<GaussianBasisSet>& bs, const Ref<MessageGrp>& g,
    00051                 signed char *pm) :
    00052       GBuild<T>(t), grp_(g), tbi_(tbi), integral_(ints), gbs_(bs), pmax(pm) {}
    00053     ~LocalLBGBuild() {}
    00054
    00055     void build_gmat(double accuracy) {
    00056       tim_enter("ao_gmat");
    00057       tim_set_default("quartet");
    00058
    00059       double tnint=0;
    00060       int tol = (int) (log(accuracy)/log(2.0));
    00061       int me=grp_->me();
    00062       int nproc = grp_->n();
    00063
    00064       Ref<PetiteList> rpl = integral_->petite_list();
    00065
    00066       // grab references for speed
    00067       GaussianBasisSet& gbs = *gbs_.pointer();
    00068       PetiteList& pl = *rpl.pointer();
    00069       TwoBodyInt& tbi = *tbi_.pointer();
    00070
    00071       tbi.set_redundant(0);
    00072       const double *intbuf = tbi.buffer();
    00073
    00074       int inds[4];
    00075
    00076       // node zero passes out indices
    00077       if (me==0) {
    00078         int i;
    00079         for (i=0; i < gbs.nshell(); i++) {
    00080           if (!pl.in_p1(i))
    00081             continue;
    00082
    00083           inds[0]=i;
    00084
    00085           for (int j=0; j <= i; j++) {
    00086             int oij = i_offset(i)+j;
    00087
    00088             if (!pl.in_p2(oij))
    00089               continue;
    00090
    00091             inds[1]=j;
    00092
    00093             tim_enter_default();
    00094             int from;
    00095             grp_->recvt(2323, &from, 1);
    00096             grp_->sendt(from, 3232, inds, 4);
    00097             tim_exit_default();
    00098           }
    00099         }
    00100
    00101         // now clean up
    00102         inds[0] = inds[1] = inds[2] = inds[3] = -1;
    00103
    00104         for (i=1; i < nproc; i++) {
    00105           int from;
    00106           grp_->recvt(2323, &from, 1);
    00107           grp_->sendt(from, 3232, inds, 4);
    00108         }
    00109       }
    00110
    00111       // all other nodes do the work
    00112       else {
    00113         do {
    00114           grp_->sendt(0, 2323, &me, 1);
    00115           grp_->recvt(3232, inds, 4);
    00116
    00117           int i=inds[0];
    00118           int j=inds[1];
    00119
    00120           if (i < 0)
    00121             break;
    00122
    00123           int fi=gbs.shell_to_function(i);
    00124           int ni=gbs(i).nfunction();
    00125
    00126           int oij = i_offset(i)+j;
    00127
    00128           int fj=gbs.shell_to_function(j);
    00129           int nj=gbs(j).nfunction();
    00130           int pmaxij = pmax[oij];
    00131
    00132           for (int k=0; k <= i; k++) {
    00133           int fk=gbs.shell_to_function(k);
    00134           int nk=gbs(k).nfunction();
    00135
    00136           int pmaxijk=pmaxij, ptmp;
    00137           if ((ptmp=pmax[i_offset(i)+k]-2) > pmaxijk) pmaxijk=ptmp;
    00138           if ((ptmp=pmax[ij_offset(j,k)]-2) > pmaxijk) pmaxijk=ptmp;
    00139
    00140           int okl = i_offset(k);
    00141           for (int l=0; l <= (k==i?j:k); l++,okl++) {
    00142
    00143             int pmaxijkl = pmaxijk;
    00144             if ((ptmp=pmax[okl]) > pmaxijkl) pmaxijkl=ptmp;
    00145             if ((ptmp=pmax[i_offset(i)+l]-2) > pmaxijkl) pmaxijkl=ptmp;
    00146             if ((ptmp=pmax[ij_offset(j,l)]-2) > pmaxijkl) pmaxijkl=ptmp;
    00147
    00148
    00149             if (tbi.log2_shell_bound(i,j,k,l)+pmaxijkl < tol)
    00150               continue;
    00151
    00152             int qijkl = pl.in_p4(oij,okl,i,j,k,l);
    00153             if (!qijkl)
    00154               continue;
    00155
    00156             tim_enter_default();
    00157             tbi.compute_shell(i,j,k,l);
    00158             tim_exit_default();
    00159
    00160             int e12 = (i==j);
    00161             int e34 = (k==l);
    00162             int e13e24 = (i==k) && (j==l);
    00163             int e_any = e12||e34||e13e24;
    00164
    00165             int fl=gbs.shell_to_function(l);
    00166             int nl=gbs(l).nfunction();
    00167
    00168             int ii,jj,kk,ll;
    00169             int I,J,K,L;
    00170             int index=0;
    00171
    00172             for (I=0, ii=fi; I < ni; I++, ii++) {
    00173               for (J=0, jj=fj; J <= (e12 ? I : nj-1); J++, jj++) {
    00174                 for (K=0, kk=fk; K <= (e13e24 ? I : nk-1); K++, kk++) {
    00175
    00176                   int lend = (e34 ? ((e13e24)&&(K==I) ? J : K)
    00177                               : ((e13e24)&&(K==I)) ? J : nl-1);
    00178                   for (L=0, ll=fl; L <= lend; L++, ll++, index++) {
    00179                     double pki_int = intbuf[index];
    00180
    00181                     if ((pki_int>0?pki_int:-pki_int) < 1.0e-15)
    00182                       continue;
    00183
    00184                     if (qijkl > 1)
    00185                       pki_int *= qijkl;
    00186
    00187       if (e_any) {
    00188         int ij,kl;
    00189         double val;
    00190
    00191         if (jj == kk) {
    00192           /*
    00193            * if i=j=k or j=k=l, then this integral contributes
    00194            * to J, K1, and K2 of G(ij), so
    00195            * pkval = (ijkl) - 0.25 * ((ikjl)-(ilkj))
    00196            *       = 0.5 * (ijkl)
    00197            */
    00198           if (ii == jj || kk == ll) {
    00199             ij = i_offset(ii)+jj;
    00200             kl = i_offset(kk)+ll;
    00201             val = (ij==kl) ? 0.5*pki_int : pki_int;
    00202
    00203             contribution.cont5(ij,kl,val);
    00204
    00205           } else {
    00206             /*
    00207              * if j=k, then this integral contributes
    00208              * to J and K1 of G(ij)
    00209              *
    00210              * pkval = (ijkl) - 0.25 * (ikjl)
    00211              *       = 0.75 * (ijkl)
    00212              */
    00213             ij = i_offset(ii)+jj;
    00214             kl = i_offset(kk)+ll;
    00215             val = (ij==kl) ? 0.5*pki_int : pki_int;
    00216
    00217             contribution.cont4(ij,kl,val);
    00218
    00219             /*
    00220              * this integral also contributes to K1 and K2 of
    00221              * G(il)
    00222              *
    00223              * pkval = -0.25 * ((ijkl)+(ikjl))
    00224              *       = -0.5 * (ijkl)
    00225              */
    00226             ij = ij_offset(ii,ll);
    00227             kl = ij_offset(kk,jj);
    00228             val = (ij==kl) ? 0.5*pki_int : pki_int;
    00229
    00230             contribution.cont3(ij,kl,val);
    00231           }
    00232         } else if (ii == kk || jj == ll) {
    00233           /*
    00234            * if i=k or j=l, then this integral contributes
    00235            * to J and K2 of G(ij)
    00236            *
    00237            * pkval = (ijkl) - 0.25 * (ilkj)
    00238            *       = 0.75 * (ijkl)
    00239            */
    00240           ij = i_offset(ii)+jj;
    00241           kl = i_offset(kk)+ll;
    00242           val = (ij==kl) ? 0.5*pki_int : pki_int;
    00243
    00244           contribution.cont4(ij,kl,val);
    00245
    00246           /*
    00247            * this integral also contributes to K1 and K2 of
    00248            * G(ik)
    00249            *
    00250            * pkval = -0.25 * ((ijkl)+(ilkj))
    00251            *       = -0.5 * (ijkl)
    00252            */
    00253           ij = ij_offset(ii,kk);
    00254           kl = ij_offset(jj,ll);
    00255           val = (ij==kl) ? 0.5*pki_int : pki_int;
    00256
    00257           contribution.cont3(ij,kl,val);
    00258
    00259         } else {
    00260           /*
    00261            * This integral contributes to J of G(ij)
    00262            *
    00263            * pkval = (ijkl)
    00264            */
    00265           ij = i_offset(ii)+jj;
    00266           kl = i_offset(kk)+ll;
    00267           val = (ij==kl) ? 0.5*pki_int : pki_int;
    00268
    00269           contribution.cont1(ij,kl,val);
    00270
    00271           /*
    00272            * and to K1 of G(ik)
    00273            *
    00274            * pkval = -0.25 * (ijkl)
    00275            */
    00276           ij = ij_offset(ii,kk);
    00277           kl = ij_offset(jj,ll);
    00278           val = (ij==kl) ? 0.5*pki_int : pki_int;
    00279
    00280           contribution.cont2(ij,kl,val);
    00281
    00282           if ((ii != jj) && (kk != ll)) {
    00283             /*
    00284              * if i!=j and k!=l, then this integral also
    00285              * contributes to K2 of G(il)
    00286              *
    00287              * pkval = -0.25 * (ijkl)
    00288              *
    00289              * note: if we get here, then ik can't equal jl,
    00290              * so pkval wasn't multiplied by 0.5 above.
    00291              */
    00292             ij = ij_offset(ii,ll);
    00293             kl = ij_offset(kk,jj);
    00294
    00295             contribution.cont2(ij,kl,val);
    00296           }
    00297         }
    00298       } else { // !e_any
    00299         if (jj == kk) {
    00300           /*
    00301            * if j=k, then this integral contributes
    00302            * to J and K1 of G(ij)
    00303            *
    00304            * pkval = (ijkl) - 0.25 * (ikjl)
    00305            *       = 0.75 * (ijkl)
    00306            */
    00307           contribution.cont4(i_offset(ii)+jj,i_offset(kk)+ll,pki_int);
    00308
    00309           /*
    00310            * this integral also contributes to K1 and K2 of
    00311            * G(il)
    00312            *
    00313            * pkval = -0.25 * ((ijkl)+(ikjl))
    00314            *       = -0.5 * (ijkl)
    00315            */
    00316           contribution.cont3(ij_offset(ii,ll),ij_offset(kk,jj),pki_int);
    00317
    00318         } else if (ii == kk || jj == ll) {
    00319           /*
    00320            * if i=k or j=l, then this integral contributes
    00321            * to J and K2 of G(ij)
    00322            *
    00323            * pkval = (ijkl) - 0.25 * (ilkj)
    00324            *       = 0.75 * (ijkl)
    00325            */
    00326           contribution.cont4(i_offset(ii)+jj,i_offset(kk)+ll,pki_int);
    00327
    00328           /*
    00329            * this integral also contributes to K1 and K2 of
    00330            * G(ik)
    00331            *
    00332            * pkval = -0.25 * ((ijkl)+(ilkj))
    00333            *       = -0.5 * (ijkl)
    00334            */
    00335           contribution.cont3(ij_offset(ii,kk),ij_offset(jj,ll),pki_int);
    00336
    00337         } else {
    00338           /*
    00339            * This integral contributes to J of G(ij)
    00340            *
    00341            * pkval = (ijkl)
    00342            */
    00343           contribution.cont1(i_offset(ii)+jj,i_offset(kk)+ll,pki_int);
    00344
    00345           /*
    00346            * and to K1 of G(ik)
    00347            *
    00348            * pkval = -0.25 * (ijkl)
    00349            */
    00350           contribution.cont2(ij_offset(ii,kk),ij_offset(jj,ll),pki_int);
    00351
    00352           /*
    00353            * and to K2 of G(il)
    00354            *
    00355            * pkval = -0.25 * (ijkl)
    00356            */
    00357           contribution.cont2(ij_offset(ii,ll),ij_offset(kk,jj),pki_int);
    00358         }
    00359       }
    00360                   }
    00361                 }
    00362               }
    00363             }
    00364
    00365             tnint += (double) ni*nj*nk*nl;
    00366           }
    00367           }
    00368         } while (inds[0] > -1);
    00369       }
    00370
    00371       grp_->sum(&tnint, 1, 0, 0);
    00372       ExEnv::out0() << indent << scprintf("%20.0f integrals\n", tnint);
    00373
    00374       tim_exit("ao_gmat");
    00375     }
    00376
    00377 };
    00378
    00379 }
    00380
    00381 #endif
    00382
    00383 // Local Variables:
    00384 // mode: c++
    00385 // c-file-style: "ETS"
    00386 // End: