obint.cc

大型并行量子化学软件；支持密度泛函（DFT）。可以进行各种量子化学计算。支持CHARMM并行计算。非常具有应用价值。

      
      if (ifn < b->istart || ifn >= b->iend ||
          jfn < b->jstart || jfn >= b->jend)
        continue;
      
      int data_index = (ifn - b->istart)*njdata + jfn - b->jstart;
      b->data[data_index] += spi.val();
    }
  }
}

void
OneBodyIntOp::process_spec_ltri(SCMatrixLTriBlock* b)
{
  Ref<GaussianBasisSet> bs1 = iter->one_body_int()->basis1();

  // convert basis function indices into shell indices
  int fnstart = b->start;
  int fnend = b->end;
  int shstart = bs1->function_to_shell(fnstart);
  int shend = (fnend?bs1->function_to_shell(fnend - 1) + 1 : 0);

  iter->set_redundant(1);

  // loop over all needed shells
  for (iter->start(shstart,shstart,shend,shend); iter->ready(); iter->next()) {
    ShellPairIter& spi = iter->current_pair();

    // compute a set of shell integrals
    for (spi.start(); spi.ready(); spi.next()) {
      int ifn = spi.i();
      int jfn = spi.j();
      
      if (ifn < fnstart || ifn >= fnend)
        continue;
      
      int ioff = ifn-fnstart;
      int joff = jfn-fnstart;
      
      int data_index = i_offset(ioff)+joff;
      
      b->data[data_index] += spi.val();
    }
  }
}

void
OneBodyIntOp::process_spec_rectsub(SCMatrixRectSubBlock* b)
{
  Ref<GaussianBasisSet> bs1 = iter->one_body_int()->basis1();
  Ref<GaussianBasisSet> bs2 = iter->one_body_int()->basis2();
  
  // convert basis function indices into shell indices
  int istart = b->istart;
  int jstart = b->jstart;
  int iend = b->iend;
  int jend = b->jend;

  int ishstart = bs1->function_to_shell(istart);
  int jshstart = bs2->function_to_shell(jstart);

  int ishend = (iend ? bs1->function_to_shell(iend-1) + 1 : 0);
  int jshend = (jend ? bs2->function_to_shell(jend-1) + 1 : 0);

  int njdata = b->istride;

  iter->set_redundant(0);

  for (iter->start(ishstart,jshstart,ishend,jshend);
       iter->ready(); iter->next()) {
    ShellPairIter& spi = iter->current_pair();

    for (spi.start(); spi.ready(); spi.next()) {
      int ifn = spi.i();
      int jfn = spi.j();
      
      if (ifn < istart || ifn >= iend || jfn < jstart || jfn >= jend)
        continue;
      
      int data_index = ifn*njdata + jfn;
      b->data[data_index] += spi.val();
    }
  }
}

void
OneBodyIntOp::process_spec_ltrisub(SCMatrixLTriSubBlock* b)
{
  Ref<GaussianBasisSet> bs1 = iter->one_body_int()->basis1();

  // convert basis function indices into shell indices
  int istart = b->istart;
  int iend = b->iend;

  int jstart = b->jstart;
  int jend = b->jend;

  int ishstart = bs1->function_to_shell(istart);
  int jshstart = bs1->function_to_shell(jstart);

  int ishend = (iend ? bs1->function_to_shell(iend-1) + 1 : 0);
  int jshend = (jend ? bs1->function_to_shell(jend-1) + 1 : 0);

  iter->set_redundant(1);

  // loop over all needed shells
  for (iter->start(ishstart,jshstart,ishend,jshend);
       iter->ready(); iter->next()) {
    ShellPairIter& spi = iter->current_pair();

    // compute a set of shell integrals
    for (spi.start(); spi.ready(); spi.next()) {
      int ifn = spi.i();
      int jfn = spi.j();
      
      if (ifn < istart || ifn >= iend || jfn < jstart || jfn >= jend)
        continue;
      
      int data_index = i_offset(ifn)+jfn;
      b->data[data_index] += spi.val();
    }
  }
}

int
OneBodyIntOp::has_side_effects()
{
  return 1;
}

///////////////////////////////////////////////////////////////////////

OneBody3IntOp::OneBody3IntOp(const Ref<OneBodyInt>& it)
{
  iter = new OneBodyIntIter(it);
}

OneBody3IntOp::OneBody3IntOp(const Ref<OneBodyIntIter>& it) :
  iter(it)
{
}

OneBody3IntOp::~OneBody3IntOp()
{
}

void
OneBody3IntOp::process(SCMatrixBlockIter&,
                       SCMatrixBlockIter&,
                       SCMatrixBlockIter&)
{
  ExEnv::err0() << indent
       << "OneBody3IntOp::process(SCMatrixBlockIter&): "
       << "cannot handle generic case\n";
  abort();
}

void
OneBody3IntOp::process_spec_rect(SCMatrixRectBlock* a,
                                 SCMatrixRectBlock* b,
                                 SCMatrixRectBlock* c)
{
  Ref<GaussianBasisSet> bs1 = iter->one_body_int()->basis1();
  Ref<GaussianBasisSet> bs2 = iter->one_body_int()->basis2();

  // convert basis function indices into shell indices
  int ishstart = bs1->function_to_shell(b->istart);
  int jshstart = bs2->function_to_shell(b->jstart);

  int ishend = bs1->function_to_shell(b->iend);
  int jshend = bs2->function_to_shell(b->jend);

  iter->set_redundant(0);

  for (iter->start(ishstart,jshstart,ishend,jshend);
       iter->ready(); iter->next()) {
    // compute a set of shell integrals
    ShellPairIter& spi = iter->current_pair();

    for (spi.start(); spi.ready(); spi.next()) {
      int ifn = spi.i();
      int jfn = spi.j();
        
      if (ifn < b->istart || ifn >= b->iend ||
          jfn < b->jstart || jfn >= b->jend)
        continue;
      
#if 0
      for (int i=0; i<nish; i++,ifn++) {
          if (ifn < b->istart || ifn >= b->iend) {
              tmp += njsh * 3;
            }
          else {
              int jfn = jfnsave;
              int data_index = (ifn - b->istart)*njdata + jfn - b->jstart;
              for (int j=0; j<njsh; j++,jfn++) {
                  if (jfn >= b->jstart && jfn < b->jend) {
                      a->data[data_index] += tmp[0] * scale;
                      b->data[data_index] += tmp[1] * scale;
                      c->data[data_index] += tmp[2] * scale;
                      data_index++;
                    }
                  tmp += 3;
                }
            }
        }
#endif
    }
  }
}

void
OneBody3IntOp::process_spec_ltri(SCMatrixLTriBlock* a,
                                 SCMatrixLTriBlock* b,
                                 SCMatrixLTriBlock* c)
{
#if 0
  Ref<GaussianBasisSet> bs1 = iter->one_body_int()->basis1();

  // convert basis function indices into shell indices
  int fnstart = b->start;
  int fnend = b->end;
  int shstart = bs1->function_to_shell(fnstart);
  int shend = (fnend?bs1->function_to_shell(fnend - 1) + 1 : 0);

  // loop over all needed shells
  iter->reset(shstart, shend, 0, 0);
  for (iter->start_ltri(); iter->ready_ltri(); iter->next_ltri()) {
      int ish=iter->ishell();
      int jsh=iter->jshell();

      int nish = bs1->operator[](ish).nfunction();
      int njsh = bs1->operator[](jsh).nfunction();

      double scale = iter->scale();
    
      // compute a set of shell integrals
      compute_shell(ish,jsh,buffer_);

      // take the integrals from buffer and put them into the LTri block
      double*tmp = buffer_;
      int ifn = bs1->shell_to_function(ish);
      int jfnsave = bs1->shell_to_function(jsh);
      for (int i=0; i<nish; i++,ifn++) {
          // skip over basis functions that are not needed
          if (ifn < fnstart || ifn >= fnend) {
              tmp += njsh * 3;
            }
          else {
              int jfn = jfnsave;
              int irelfn = ifn - fnstart;
              int data_index = ((irelfn+1)*irelfn>>1) + jfn - fnstart;
              for (int j=0; j<njsh; j++,jfn++) {
                  // skip over basis functions that are not needed
                  if (jfn <= ifn && jfn >= fnstart) {
                      a->data[data_index] += tmp[0] * scale;
                      b->data[data_index] += tmp[1] * scale;
                      c->data[data_index] += tmp[2] * scale;
                      data_index++;
                    }
                  tmp += 3;
                }
            }
        }
    }
#endif
}

int
OneBody3IntOp::has_side_effects()
{
  return 1;
}

int
OneBody3IntOp::has_side_effects_in_arg1()
{
  return 1;
}

int
OneBody3IntOp::has_side_effects_in_arg2()
{
  return 1;
}

///////////////////////////////////////////////////////////////////////

OneBodyDerivInt::OneBodyDerivInt(Integral *integral,
                                 const Ref<GaussianBasisSet>&b) :
  integral_(integral),
  bs1(b), bs2(b)
{
  // allocate a buffer
  int biggest_shell = b->max_nfunction_in_shell();
  biggest_shell *= biggest_shell * 3;

  if (biggest_shell) {
    buffer_ = new double[biggest_shell];
  } else {
    buffer_ = 0;
  }
}

OneBodyDerivInt::OneBodyDerivInt(Integral *integral,
                                 const Ref<GaussianBasisSet>&b1,
                                 const Ref<GaussianBasisSet>&b2) :
  integral_(integral),
  bs1(b1), bs2(b2)
{
  buffer_ = 0;
}

OneBodyDerivInt::~OneBodyDerivInt()
{
}

int
OneBodyDerivInt::nbasis() const
{
  return bs1->nbasis();
}

int
OneBodyDerivInt::nbasis1() const
{
  return bs1->nbasis();
}

int
OneBodyDerivInt::nbasis2() const
{
  return bs2->nbasis();
}

int
OneBodyDerivInt::nshell() const
{
  return bs1->nshell();
}

int
OneBodyDerivInt::nshell1() const
{
  return bs1->nshell();
}

int
OneBodyDerivInt::nshell2() const
{
  return bs2->nshell();
}

Ref<GaussianBasisSet>
OneBodyDerivInt::basis()
{
  return bs1;
}

Ref<GaussianBasisSet>
OneBodyDerivInt::basis1()
{
  return bs1;
}

Ref<GaussianBasisSet>
OneBodyDerivInt::basis2()
{
  return bs2;
}

const double *
OneBodyDerivInt::buffer() const
{
  return buffer_;
}

/////////////////////////////////////////////////////////////////////////////

// Local Variables:
// mode: c++
// c-file-style: "ETS"
// End: