nonbonded.c

最著名最快的分子模拟软件

/*
 * $Id: nonbonded.c,v 1.14.2.5 2008/02/29 07:02:44 spoel Exp $
 * 
 *                This source code is part of
 * 
 *                 G   R   O   M   A   C   S
 * 
 *          GROningen MAchine for Chemical Simulations
 * 
 *                        VERSION 3.3.3
 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
 * Copyright (c) 2001-2008, The GROMACS development team,
 * check out http://www.gromacs.org for more information.

 * This program is free software; you can redistribute it and/or
 * modify it under the terms of the GNU General Public License
 * as published by the Free Software Foundation; either version 2
 * of the License, or (at your option) any later version.
 * 
 * If you want to redistribute modifications, please consider that
 * scientific software is very special. Version control is crucial -
 * bugs must be traceable. We will be happy to consider code for
 * inclusion in the official distribution, but derived work must not
 * be called official GROMACS. Details are found in the README & COPYING
 * files - if they are missing, get the official version at www.gromacs.org.
 * 
 * To help us fund GROMACS development, we humbly ask that you cite
 * the papers on the package - you can find them in the top README file.
 * 
 * For more info, check our website at http://www.gromacs.org
 * 
 * And Hey:
 * Groningen Machine for Chemical Simulation
 */
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif


#include <gmx_thread.h>


#include <stdio.h>
#include <stdlib.h>
#include "typedefs.h"
#include "txtdump.h"
#include "smalloc.h"
#include "ns.h"
#include "vec.h"
#include "maths.h"
#include "macros.h"
#include "force.h"
#include "names.h"
#include "main.h"
#include "xvgr.h"
#include "gmx_fatal.h"
#include "physics.h"
#include "force.h"
#include "bondf.h"
#include "nrnb.h"
#include "smalloc.h"

#include "nb_kernel/nb_kernel.h"
#include "nb_kernel/nb_kernel330.h"
#include "nb_free_energy.h"


#ifdef GMX_PPC_ALTIVEC   
#include "nb_kernel_ppc_altivec/nb_kernel_ppc_altivec.h"
#endif

#ifdef GMX_IA32_3DNOW  
#include "nb_kernel_ia32_3dnow/nb_kernel_ia32_3dnow.h"
#endif

#ifdef GMX_IA32_SSE   
#include "nb_kernel_ia32_sse/nb_kernel_ia32_sse.h"
#endif

#ifdef GMX_IA32_SSE2    
#include "nb_kernel_ia32_sse2/nb_kernel_ia32_sse2.h"
#endif

#ifdef GMX_X86_64_SSE    
#include "nb_kernel_x86_64_sse/nb_kernel_x86_64_sse.h"
#endif

#ifdef GMX_X86_64_SSE2    
#include "nb_kernel_x86_64_sse2/nb_kernel_x86_64_sse2.h"
#endif

#if (defined GMX_IA64_ASM && defined GMX_DOUBLE) 
#include "nb_kernel_ia64_double/nb_kernel_ia64_double.h"
#endif

#if (defined GMX_IA64_ASM && !defined GMX_DOUBLE)
#include "nb_kernel_ia64_single/nb_kernel_ia64_single.h"
#endif




enum { TABLE_NONE, TABLE_COMBINED, TABLE_COUL, TABLE_VDW, TABLE_NR };


/* Table version for each kernel.
 */
static const int 
nb_kernel_table[eNR_NBKERNEL_NR] = 
{
  TABLE_NONE,     /* kernel010 */
  TABLE_NONE,     /* kernel020 */
  TABLE_VDW,      /* kernel030 */
  TABLE_NONE,     /* kernel100 */
  TABLE_NONE,     /* kernel101 */
  TABLE_NONE,     /* kernel102 */
  TABLE_NONE,     /* kernel103 */
  TABLE_NONE,     /* kernel104 */
  TABLE_NONE,     /* kernel110 */
  TABLE_NONE,     /* kernel111 */
  TABLE_NONE,     /* kernel112 */
  TABLE_NONE,     /* kernel113 */
  TABLE_NONE,     /* kernel114 */
  TABLE_NONE,     /* kernel120 */
  TABLE_NONE,     /* kernel121 */
  TABLE_NONE,     /* kernel122 */
  TABLE_NONE,     /* kernel123 */
  TABLE_NONE,     /* kernel124 */
  TABLE_VDW,      /* kernel130 */
  TABLE_VDW,      /* kernel131 */
  TABLE_VDW,      /* kernel132 */
  TABLE_VDW,      /* kernel133 */
  TABLE_VDW,      /* kernel134 */
  TABLE_NONE,     /* kernel200 */
  TABLE_NONE,     /* kernel201 */
  TABLE_NONE,     /* kernel202 */
  TABLE_NONE,     /* kernel203 */
  TABLE_NONE,     /* kernel204 */
  TABLE_NONE,     /* kernel210 */
  TABLE_NONE,     /* kernel211 */
  TABLE_NONE,     /* kernel212 */
  TABLE_NONE,     /* kernel213 */
  TABLE_NONE,     /* kernel214 */
  TABLE_NONE,     /* kernel220 */
  TABLE_NONE,     /* kernel221 */
  TABLE_NONE,     /* kernel222 */
  TABLE_NONE,     /* kernel223 */
  TABLE_NONE,     /* kernel224 */
  TABLE_VDW,      /* kernel230 */
  TABLE_VDW,      /* kernel231 */
  TABLE_VDW,      /* kernel232 */
  TABLE_VDW,      /* kernel233 */
  TABLE_VDW,      /* kernel234 */
  TABLE_COUL,     /* kernel300 */
  TABLE_COUL,     /* kernel301 */
  TABLE_COUL,     /* kernel302 */
  TABLE_COUL,     /* kernel303 */
  TABLE_COUL,     /* kernel304 */
  TABLE_COUL,     /* kernel310 */
  TABLE_COUL,     /* kernel311 */
  TABLE_COUL,     /* kernel312 */
  TABLE_COUL,     /* kernel313 */
  TABLE_COUL,     /* kernel314 */
  TABLE_COUL,     /* kernel320 */
  TABLE_COUL,     /* kernel321 */
  TABLE_COUL,     /* kernel322 */
  TABLE_COUL,     /* kernel323 */
  TABLE_COUL,     /* kernel324 */
  TABLE_COMBINED, /* kernel330 */
  TABLE_COMBINED, /* kernel331 */
  TABLE_COMBINED, /* kernel332 */
  TABLE_COMBINED, /* kernel333 */
  TABLE_COMBINED, /* kernel334 */
  TABLE_NONE,     /* kernel400 */
  TABLE_NONE,     /* kernel410 */
  TABLE_VDW       /* kernel430 */
};



static nb_kernel_t **
nb_kernel_list = NULL;


static void
setup_kernels(FILE *log)
{
    int i;
    
    snew(nb_kernel_list,eNR_NBKERNEL_NR);
    
    /* Note that later calls overwrite earlier, so the preferred (fastest)
     * version should be at the end. For instance, we call SSE after 3DNow.
     */
    
    for(i=0;i<eNR_NBKERNEL_NR;i++)
        nb_kernel_list[i] = NULL;
    
	if(log)
	    fprintf(log,"Configuring nonbonded kernels...\n");
	
    nb_kernel_setup(log,nb_kernel_list);

    if(getenv("NOASSEMBLYLOOPS") != NULL)
    {
        if(log)
            fprintf(log,
                    "Found environment variable NOASSEMBLYLOOPS.\n"
                    "Disabling all SSE/SSE2/3DNow/Altivec/ia64 asm support.\n\n");
        return;
    }
    
#ifdef GMX_DOUBLE

    /* Double precision */    

#ifdef GMX_IA32_SSE2    
    nb_kernel_setup_ia32_sse2(log,nb_kernel_list);
#elif defined GMX_X86_64_SSE2 
    nb_kernel_setup_x86_64_sse2(log,nb_kernel_list);
#elif defined GMX_IA64_ASM
    nb_kernel_setup_ia64_double(log,nb_kernel_list);
#endif
    
#else
    /* Single */

    /* Older Athlons only support 3DNow, so check that first, and 
     * then newer AthlonXP/Opteron CPUs with SSE support will override
     * it in the SSE check.
     */   
#ifdef GMX_IA32_3DNOW
    nb_kernel_setup_ia32_3dnow(log,nb_kernel_list);
#endif
#ifdef GMX_IA32_SSE    
    nb_kernel_setup_ia32_sse(log,nb_kernel_list);
#elif defined GMX_X86_64_SSE   
    nb_kernel_setup_x86_64_sse(log,nb_kernel_list);
#elif defined GMX_PPC_ALTIVEC
    nb_kernel_setup_ppc_altivec(log,nb_kernel_list);
#elif defined GMX_IA64_ASM
    nb_kernel_setup_ia64_single(log,nb_kernel_list);
#endif
    
#endif /* precision */

	if(log)
	    fprintf(log,"\n\n");
}


void do_nonbonded(FILE *fplog,t_commrec *cr,t_forcerec *fr,
                  rvec x[],rvec f[],t_mdatoms *mdatoms,
                  real egnb[],real egcoul[],rvec box_size,
                  t_nrnb *nrnb,real lambda,real *dvdlambda,
                  bool bLR,int nls,int eNL,bool bDoForces)
{
  t_nblist *      nlist;
  real *          fshift;
  int             n,n0,n1,i,i0,i1,nrnb_ind,sz;
  t_nblists       *nblists;
  bool            bWater;
  nb_kernel_t *   kernelptr;
  FILE *          fp;
  int             wateratoms;
  int             nthreads = 1;
  int             tabletype;
  int             outeriter,inneriter;
  real *          tabledata = NULL;

  if(fr->solvent_opt == esolSPC)
    {
    wateratoms = 3;
    }
  else if(fr->solvent_opt == esolTIP4P)
    {
    wateratoms = 4;
    }
  else
    {
    wateratoms = 1;
    }
  
    if (eNL >= 0) 
    {
    i0 = eNL;
    i1 = i0+1;
    }
    else
    {
    i0 = 0;
    i1 = eNL_NR;
  }
  
  if (nls >= 0) {
    n0 = nls;
    n1 = nls+1;
  } else {
    n0 = 0;
    n1 = fr->nnblists;
  }
  
  if(nb_kernel_list == NULL)
    {
    setup_kernels(fplog);    
    }
  
  if (bLR)
    {
    fshift = fr->fshift_twin[0];
    }
  else
    {
    fshift = fr->fshift[0];
    }
  
  for(n=n0; (n<n1); n++) {
    nblists = &fr->nblists[n];
    for(i=i0; (i<i1); i++) {
      outeriter = inneriter = 0;
      
      if (bLR) 
	nlist = &(nblists->nlist_lr[i]);
      else
	nlist = &(nblists->nlist_sr[i]);
      
        if (nlist->nri > 0) 
        {
	nrnb_ind = nlist->il_code;

            if(nrnb_ind==eNR_NBKERNEL_FREE_ENERGY)
            {
                /* generic free energy, use combined table */
  	        tabledata = nblists->tab.tab;
            }
            else
            {
	        tabletype = nb_kernel_table[nrnb_ind];
	      
                /* normal kernels, not free energy */
                if (!bDoForces)
                    nrnb_ind += eNR_NBKERNEL_NR/2;
                
                if(tabletype == TABLE_COMBINED)
                {
                    tabledata = nblists->tab.tab;
                }
                else if(tabletype == TABLE_COUL)
                {
                    tabledata = nblists->coultab;
                }
                else if(tabletype == TABLE_VDW)
                {
                    tabledata = nblists->vdwtab;
                }
                else
                {
                    tabledata = NULL;
                }
            }
            
            nlist->count = 0;
            
            if(nlist->free_energy)
            {
                if(nlist->ivdw==2)
                {
                    gmx_fatal(FARGS,"Cannot do free energy Buckingham interactions.");
                }