• 设为首页
  • 点击收藏
  • 手机版
    手机扫一扫访问
    迪恩网络手机版
  • 关注官方公众号
    微信扫一扫关注
    公众号

C++ rvec_inc函数代码示例

原作者: [db:作者] 来自: [db:来源] 收藏 邀请

本文整理汇总了C++中rvec_inc函数的典型用法代码示例。如果您正苦于以下问题:C++ rvec_inc函数的具体用法?C++ rvec_inc怎么用?C++ rvec_inc使用的例子?那么恭喜您, 这里精选的函数代码示例或许可以为您提供帮助。



在下文中一共展示了rvec_inc函数的20个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于我们的系统推荐出更棒的C++代码示例。

示例1: calc_com2

/* calculates com of all atoms in x[], *index has their index numbers
   to get the masses from atom[] */
real calc_com2(rvec x[],int gnx,atom_id *index,t_mdatoms *md,rvec com,
	       matrix box)
{
  int  i,ii,m;
  real m0,tm;

  clear_rvec(com);
  tm=0;
  for(i=0; (i<gnx); i++) {
    ii=index[i];
    m0=md->massT[ii];
    tm+=m0;
    for(m=0; (m<DIM); m++)
      com[m]+=m0*x[i][m];
  }
  svmul(1/tm,com,com);
  for(m=DIM-1; m>=0; m--) {
    /* next two lines used to be commented out */
    if (com[m] < 0        ) rvec_inc(com,box[m]);
    if (com[m] > box[m][m]) rvec_dec(com,box[m]); 
  }
  return tm;
}
开发者ID:Chadi-akel,项目名称:cere,代码行数:25,代码来源:pullutil.c


示例2: orient_princ

void orient_princ(t_atoms *atoms, int isize, atom_id *index,
                  int natoms, rvec x[], rvec *v, rvec d)
{
    int     i, m;
    rvec    xcm, prcomp;
    matrix  trans;

    calc_xcm(x, isize, index, atoms->atom, xcm, FALSE);
    for (i = 0; i < natoms; i++)
    {
        rvec_dec(x[i], xcm);
    }
    principal_comp(isize, index, atoms->atom, x, trans, prcomp);
    if (d)
    {
        copy_rvec(prcomp, d);
    }

    /* Check whether this trans matrix mirrors the molecule */
    if (det(trans) < 0)
    {
        for (m = 0; (m < DIM); m++)
        {
            trans[ZZ][m] = -trans[ZZ][m];
        }
    }
    rotate_atoms(natoms, NULL, x, trans);
    if (v)
    {
        rotate_atoms(natoms, NULL, v, trans);
    }

    for (i = 0; i < natoms; i++)
    {
        rvec_inc(x[i], xcm);
    }
}
开发者ID:dkarkoulis,项目名称:gromacs,代码行数:37,代码来源:princ.c


示例3: gmx_pme_receive_f

void gmx_pme_receive_f(t_commrec *cr,
                       rvec f[], matrix vir_q, real *energy_q,
                       matrix vir_lj, real *energy_lj,
                       real *dvdlambda_q, real *dvdlambda_lj,
                       float *pme_cycles)
{
    int natoms, i;

#ifdef GMX_PME_DELAYED_WAIT
    /* Wait for the x request to finish */
    gmx_pme_send_coeffs_coords_wait(cr->dd);
#endif

    natoms = cr->dd->nat_home;

    if (natoms > cr->dd->pme_recv_f_alloc)
    {
        cr->dd->pme_recv_f_alloc = over_alloc_dd(natoms);
        srenew(cr->dd->pme_recv_f_buf, cr->dd->pme_recv_f_alloc);
    }

#ifdef GMX_MPI
    MPI_Recv(cr->dd->pme_recv_f_buf[0],
             natoms*sizeof(rvec), MPI_BYTE,
             cr->dd->pme_nodeid, 0, cr->mpi_comm_mysim,
             MPI_STATUS_IGNORE);
#endif

    for (i = 0; i < natoms; i++)
    {
        rvec_inc(f[i], cr->dd->pme_recv_f_buf[i]);
    }


    receive_virial_energy(cr, vir_q, energy_q, vir_lj, energy_lj, dvdlambda_q, dvdlambda_lj, pme_cycles);
}
开发者ID:alfredog,项目名称:gromacs,代码行数:36,代码来源:pme-pp.cpp


示例4: reduce_thread_forces

static void reduce_thread_forces(int n, rvec *f,
                                 tensor vir,
                                 real *Vcorr,
                                 int efpt_ind, real *dvdl,
                                 int nthreads, f_thread_t *f_t)
{
    int t, i;

    /* This reduction can run over any number of threads */
#pragma omp parallel for num_threads(gmx_omp_nthreads_get(emntBonded)) private(t) schedule(static)
    for (i = 0; i < n; i++)
    {
        for (t = 1; t < nthreads; t++)
        {
            rvec_inc(f[i], f_t[t].f[i]);
        }
    }
    for (t = 1; t < nthreads; t++)
    {
        *Vcorr += f_t[t].Vcorr;
        *dvdl  += f_t[t].dvdl[efpt_ind];
        m_add(vir, f_t[t].vir, vir);
    }
}
开发者ID:yhalcyon,项目名称:Gromacs,代码行数:24,代码来源:force.c


示例5: calc_geom

real calc_geom(int isize, atom_id *index, rvec *x, rvec geom_center, rvec min,
               rvec max, gmx_bool bDiam)
{
    real  diam2, d;
    char *grpnames;
    int   ii, i, j;

    clear_rvec(geom_center);
    diam2 = 0;
    if (isize == 0)
    {
        clear_rvec(min);
        clear_rvec(max);
    }
    else
    {
        if (index)
        {
            ii = index[0];
        }
        else
        {
            ii = 0;
        }
        for (j = 0; j < DIM; j++)
        {
            min[j] = max[j] = x[ii][j];
        }
        for (i = 0; i < isize; i++)
        {
            if (index)
            {
                ii = index[i];
            }
            else
            {
                ii = i;
            }
            rvec_inc(geom_center, x[ii]);
            for (j = 0; j < DIM; j++)
            {
                if (x[ii][j] < min[j])
                {
                    min[j] = x[ii][j];
                }
                if (x[ii][j] > max[j])
                {
                    max[j] = x[ii][j];
                }
            }
            if (bDiam)
            {
                if (index)
                {
                    for (j = i + 1; j < isize; j++)
                    {
                        d     = distance2(x[ii], x[index[j]]);
                        diam2 = max(d, diam2);
                    }
                }
                else
                {
                    for (j = i + 1; j < isize; j++)
                    {
                        d     = distance2(x[i], x[j]);
                        diam2 = max(d, diam2);
                    }
                }
            }
        }
        svmul(1.0 / isize, geom_center, geom_center);
    }

    return sqrt(diam2);
}
开发者ID:yhalcyon,项目名称:Gromacs,代码行数:75,代码来源:gmx_editconf.c


示例6: list_trn

static void list_trn(char *fn)
{
  static real mass[5] = { 15.9994, 1.008, 1.008, 0.0, 0.0 };
  int         i,j=0,m,fpread,fpwrite,nframe;
  rvec        *x,*v,*f,fmol[2],xcm[2],torque[j],dx;
  real        mmm,len;
  matrix      box;
  t_trnheader trn;
  gmx_bool        bOK;

  printf("Going to open %s\n",fn);
  fpread  = open_trn(fn,"r"); 
  fpwrite = open_tpx(NULL,"w");
  gmx_fio_setdebug(fpwrite,TRUE);
  
  mmm=mass[0]+2*mass[1];
  for(i=0; (i<5); i++) 
    mass[i] /= mmm;
  
  nframe = 0;
  while (fread_trnheader(fpread,&trn,&bOK)) {
    snew(x,trn.natoms);
    snew(v,trn.natoms);
    snew(f,trn.natoms);
    if (fread_htrn(fpread,&trn,
		   trn.box_size ? box : NULL,
		   trn.x_size   ? x : NULL,
		   trn.v_size   ? v : NULL,
		   trn.f_size   ? f : NULL)) {
		   
      if (trn.x_size && trn.f_size) {
	printf("There are %d atoms\n",trn.natoms);
	for(j=0; (j<2); j++) {
	  clear_rvec(xcm[j]);
	  clear_rvec(fmol[j]);
	  clear_rvec(torque[j]);
	  for(i=5*j; (i<5*j+5); i++) {
	    rvec_inc(fmol[j],f[i]);
	    for(m=0; (m<DIM); m++)
	      xcm[j][m] += mass[i%5]*x[i][m];
	  }
	  for(i=5*j; (i<5*j+5); i++) {
	    rvec_dec(x[i],xcm[j]);
	    cprod(x[i],f[i],dx);
	    rvec_inc(torque[j],dx);
	    rvec_inc(x[i],xcm[j]);
	  }
	}
	pr_rvecs(stdout,0,"FMOL  ",fmol,2);
	pr_rvecs(stdout,0,"TORQUE",torque,2);
	printf("Distance matrix Water1-Water2\n%5s","");
	for(j=0; (j<5); j++) 
	  printf("  %10s",nm[j]);
	printf("\n");
	for(j=0; (j<5); j++) {
	  printf("%5s",nm[j]);
	  for(i=5; (i<10); i++) {
	    rvec_sub(x[i],x[j],dx);
	    len = sqrt(iprod(dx,dx));
	    printf("  %10.7f",len);
	  }
	  printf("\n");
	}
      }
    }
    sfree(x);
    sfree(v);
    sfree(f);
    nframe++;
  }
  if (!bOK)
    fprintf(stderr,"\nWARNING: Incomplete frame header: nr %d, t=%g\n",
	    nframe,trn.t);
  close_tpx(fpwrite);
  close_trn(fpread);
}
开发者ID:Ruyk,项目名称:gromacs,代码行数:76,代码来源:anaf.c


示例7: gmx_traj


//.........这里部分代码省略.........
        if (bOT && fr.bV)
        {
            fprintf(outt, " %g", time);
            for (i = 0; i < ngroups; i++)
            {
                fprintf(outt, sffmt, temp(fr.v, mass, isize[i], index[i]));
            }
            fprintf(outt, "\n");
        }
        if (bEKT && fr.bV)
        {
            fprintf(outekt, " %g", time);
            for (i = 0; i < ngroups; i++)
            {
                fprintf(outekt, sffmt, ektrans(fr.v, mass, isize[i], index[i]));
            }
            fprintf(outekt, "\n");
        }
        if (bEKR && fr.bX && fr.bV)
        {
            fprintf(outekr, " %g", time);
            for (i = 0; i < ngroups; i++)
            {
                fprintf(outekr, sffmt, ekrot(fr.x, fr.v, mass, isize[i], index[i]));
            }
            fprintf(outekr, "\n");
        }
        if ((bCV || bCF) && fr.bX &&
            (ctime < 0 || (fr.time >= ctime*0.999999 &&
                           fr.time <= ctime*1.000001)))
        {
            for (i = 0; i < fr.natoms; i++)
            {
                rvec_inc(sumx[i], fr.x[i]);
            }
            nr_xfr++;
        }
        if (bCV && fr.bV)
        {
            for (i = 0; i < fr.natoms; i++)
            {
                rvec_inc(sumv[i], fr.v[i]);
            }
            nr_vfr++;
        }
        if (bCF && fr.bF)
        {
            for (i = 0; i < fr.natoms; i++)
            {
                rvec_inc(sumf[i], fr.f[i]);
            }
            nr_ffr++;
        }

    }
    while (read_next_frame(oenv, status, &fr));

    if (gpbc != NULL)
    {
        gmx_rmpbc_done(gpbc);
    }

    /* clean up a bit */
    close_trj(status);

    if (bOX)
开发者ID:MelroLeandro,项目名称:gromacs,代码行数:67,代码来源:gmx_traj.cpp


示例8: gmx_covar


//.........这里部分代码省略.........
  }
  
  /* Prepare reference frame */
  if (bPBC) {
    gpbc = gmx_rmpbc_init(&top.idef,ePBC,atoms->nr,box);
    gmx_rmpbc(gpbc,atoms->nr,box,xref);
  }
  if (bFit)
    reset_x(nfit,ifit,atoms->nr,NULL,xref,w_rls);

  snew(x,natoms);
  snew(xav,natoms);
  ndim=natoms*DIM;
  if (sqrt(GMX_LARGE_INT_MAX)<ndim) {
    gmx_fatal(FARGS,"Number of degrees of freedoms to large for matrix.\n");
  }
  snew(mat,ndim*ndim);

  fprintf(stderr,"Calculating the average structure ...\n");
  nframes0 = 0;
  nat=read_first_x(oenv,&status,trxfile,&t,&xread,box);
  if (nat != atoms->nr)
    fprintf(stderr,"\nWARNING: number of atoms in tpx (%d) and trajectory (%d) do not match\n",natoms,nat);
  do {
    nframes0++;
    /* calculate x: a fitted struture of the selected atoms */
    if (bPBC)
      gmx_rmpbc(gpbc,nat,box,xread);
    if (bFit) {
      reset_x(nfit,ifit,nat,NULL,xread,w_rls);
      do_fit(nat,w_rls,xref,xread);
    }
    for (i=0; i<natoms; i++)
      rvec_inc(xav[i],xread[index[i]]);
  } while (read_next_x(oenv,status,&t,nat,xread,box));
  close_trj(status);
  
  inv_nframes = 1.0/nframes0;
  for(i=0; i<natoms; i++)
    for(d=0; d<DIM; d++) {
      xav[i][d] *= inv_nframes;
      xread[index[i]][d] = xav[i][d];
    }
  write_sto_conf_indexed(opt2fn("-av",NFILE,fnm),"Average structure",
			 atoms,xread,NULL,epbcNONE,zerobox,natoms,index);
  sfree(xread);

  fprintf(stderr,"Constructing covariance matrix (%dx%d) ...\n",(int)ndim,(int)ndim);
  nframes=0;
  nat=read_first_x(oenv,&status,trxfile,&t,&xread,box);
  tstart = t;
  do {
    nframes++;
    tend = t;
    /* calculate x: a (fitted) structure of the selected atoms */
    if (bPBC)
      gmx_rmpbc(gpbc,nat,box,xread);
    if (bFit) {
      reset_x(nfit,ifit,nat,NULL,xread,w_rls);
      do_fit(nat,w_rls,xref,xread);
    }
    if (bRef)
      for (i=0; i<natoms; i++)
	rvec_sub(xread[index[i]],xref[index[i]],x[i]);
    else
      for (i=0; i<natoms; i++)
开发者ID:BradleyDickson,项目名称:ABPenabledGROMACS,代码行数:67,代码来源:gmx_covar.c


示例9: dd_move_f_specat

void dd_move_f_specat(gmx_domdec_t *dd, gmx_domdec_specat_comm_t *spac,
                      rvec *f, rvec *fshift)
{
    gmx_specatsend_t *spas;
    rvec             *vbuf;
    int               n, n0, n1, d, dim, dir, i;
    ivec              vis;
    int               is;
    gmx_bool          bPBC, bScrew;

    n = spac->at_end;
    for (d = dd->ndim-1; d >= 0; d--)
    {
        dim = dd->dim[d];
        if (dd->nc[dim] > 2)
        {
            /* Pulse the grid forward and backward */
            spas = spac->spas[d];
            n0   = spas[0].nrecv;
            n1   = spas[1].nrecv;
            n   -= n1 + n0;
            vbuf = spac->vbuf;
            /* Send and receive the coordinates */
            dd_sendrecv2_rvec(dd, d,
                              f+n+n1, n0, vbuf, spas[0].nsend,
                              f+n, n1, vbuf+spas[0].nsend, spas[1].nsend);
            for (dir = 0; dir < 2; dir++)
            {
                bPBC   = ((dir == 0 && dd->ci[dim] == 0) ||
                          (dir == 1 && dd->ci[dim] == dd->nc[dim]-1));
                bScrew = (bPBC && dd->bScrewPBC && dim == XX);

                spas = &spac->spas[d][dir];
                /* Sum the buffer into the required forces */
                if (!bPBC || (!bScrew && fshift == NULL))
                {
                    for (i = 0; i < spas->nsend; i++)
                    {
                        rvec_inc(f[spas->a[i]], *vbuf);
                        vbuf++;
                    }
                }
                else
                {
                    clear_ivec(vis);
                    vis[dim] = (dir == 0 ? 1 : -1);
                    is       = IVEC2IS(vis);
                    if (!bScrew)
                    {
                        /* Sum and add to shift forces */
                        for (i = 0; i < spas->nsend; i++)
                        {
                            rvec_inc(f[spas->a[i]], *vbuf);
                            rvec_inc(fshift[is], *vbuf);
                            vbuf++;
                        }
                    }
                    else
                    {
                        /* Rotate the forces */
                        for (i = 0; i < spas->nsend; i++)
                        {
                            f[spas->a[i]][XX] += (*vbuf)[XX];
                            f[spas->a[i]][YY] -= (*vbuf)[YY];
                            f[spas->a[i]][ZZ] -= (*vbuf)[ZZ];
                            if (fshift)
                            {
                                rvec_inc(fshift[is], *vbuf);
                            }
                            vbuf++;
                        }
                    }
                }
            }
        }
        else
        {
            /* Two cells, so we only need to communicate one way */
            spas = &spac->spas[d][0];
            n   -= spas->nrecv;
            /* Send and receive the coordinates */
            dd_sendrecv_rvec(dd, d, dddirForward,
                             f+n, spas->nrecv, spac->vbuf, spas->nsend);
            /* Sum the buffer into the required forces */
            if (dd->bScrewPBC && dim == XX &&
                (dd->ci[dim] == 0 ||
                 dd->ci[dim] == dd->nc[dim]-1))
            {
                for (i = 0; i < spas->nsend; i++)
                {
                    /* Rotate the force */
                    f[spas->a[i]][XX] += spac->vbuf[i][XX];
                    f[spas->a[i]][YY] -= spac->vbuf[i][YY];
                    f[spas->a[i]][ZZ] -= spac->vbuf[i][ZZ];
                }
            }
            else
            {
                for (i = 0; i < spas->nsend; i++)
                {
//.........这里部分代码省略.........
开发者ID:MelroLeandro,项目名称:gromacs,代码行数:101,代码来源:domdec_specatomcomm.cpp


示例10: do_listed_vdw_q


//.........这里部分代码省略.........
             */
            if(ivdw==2)
            {
                gmx_fatal(FARGS,"Cannot do free energy Buckingham interactions.");
            }
            count = 0;
            gmx_nb_free_energy_kernel(icoul,
                                      ivdw,
                                      i1,
                                      &i0,
                                      j_index,
                                      &i1,
                                      &shift_f,
                                      fr->shift_vec[0],
                                      fshift[0],
                                      &gid,
                                      x14[0],
                                      f14[0],
                                      chargeA,
                                      chargeB,
                                      eps,
                                      krf,
                                      crf,
                                      fr->ewaldcoeff,
                                      egcoul,
                                      typeA,
                                      typeB,
                                      ntype,
                                      nbfp,
                                      egnb,
                                      tabscale,
                                      tab,
                                      lambda,
                                      dvdlambda,
                                      fr->sc_alpha,
                                      fr->sc_power,
                                      fr->sc_sigma6_def,
                                      fr->sc_sigma6_min,
                                      TRUE,
                                      &outeriter,
                                      &inneriter);
        }
        else 
        { 
            /* Not perturbed - call kernel 330 */
            nb_kernel330
                ( &i1,
                  &i0,
                  j_index,
                  &i1,
                  &shift_f,
                  fr->shift_vec[0],
                  fshift[0],
                  &gid,
                  x14[0],
                  f14[0],
                  chargeA,
                  &eps,
                  &krf,
                  &crf,
                  egcoul,
                  typeA,
                  &ntype,
                  nbfp,
                  egnb,
                  &tabscale,
                  tab,
                  NULL,
                  NULL,
                  NULL,
                  NULL,
                  &nthreads,
                  &count,
                  (void *)&mtx,
                  &outeriter,
                  &inneriter,
                  NULL);                
        }
        
        /* Add the forces */
        rvec_inc(f[ai],f14[0]);
        rvec_dec(f[aj],f14[0]);

	if (pf_global->bInitialized)
	    pf_atom_add_bonded(pf_global, ai, aj, PF_INTER_NB14, f14[0]);

        if (g) 
        {
            /* Correct the shift forces using the graph */
            ivec_sub(SHIFT_IVEC(g,ai),SHIFT_IVEC(g,aj),dt);    
            shift_vir = IVEC2IS(dt);
            rvec_inc(fshift[shift_vir],f14[0]);
            rvec_dec(fshift[CENTRAL],f14[0]);
        }
        
	    /* flops: eNR_KERNEL_OUTER + eNR_KERNEL330 + 12 */
        }
    }
    return 0.0;
}
开发者ID:chenleo,项目名称:gromacs453pf,代码行数:101,代码来源:nonbonded.c


示例11: unitcell

void unitcell(rvec x[],rvec box,gmx_bool bYaw,real odist,real hdist)
{
#define cx  0.81649658
#define cy  0.47140452
#define cy2 0.94280904
#define cz  0.33333333
  
  rvec xx[24] = {
    { 0,   0,         0 }, /* O1 */
    { 0,   0,         1 }, /* H relative to Oxygen */
    { cx, cy,       -cz },
    { cx, cy,       -cz }, /* O2 */
    { 0, 0,       -1    }, /* H relative to Oxygen */
    { cx,-cy,       +cz },
    { cx, cy+cy2,     0 }, /* O3 */
    { -cx, cy,    -cz   }, /* H relative to Oxygen */
    { 0,   -cy2,    -cz },
    { 0,  2*cy+cy2, -cz }, /* O4 */
    {-cx,-cy,       +cz }, /* H relative to Oxygen */
    { 0 , cy2,      +cz },
    { 0,   0,         1 }, /* O5 */
    {-cx, cy,       +cz }, /* H relative to Oxygen */
    { 0 , -cy2,     +cz },
    { cx, cy,      1+cz }, /* O6 */
    { -cx, -cy,     -cz }, /* H relative to Oxygen */
    { 0,   cy2,     -cz },
    { cx, cy+cy2,     1 }, /* O7 */
    { 0,  0,       -1   }, /* H relative to Oxygen */
    { cx, cy,       +cz },
    { 0,  2*cy+cy2,1+cz }, /* O8 */
    { 0,   0,         1 }, /* H relative to Oxygen */
    { cx,   -cy,    -cz }
  };
  int  i,iin,iout,j,m;
  rvec tmp,t2,dip;
  
  clear_rvec(dip);
  for(i=0; (i<8); i++) {
    iin = 3*i;
    if (bYaw)
      iout = 5*i;
    else
      iout = iin;
    svmul(odist,xx[iin],x[iout]);
    svmul(-0.82,x[iout],t2);
    rvec_inc(dip,t2);
    for(j=1; (j<=2); j++) {
      svmul(hdist,xx[iin+j],tmp);
      rvec_add(x[iout],tmp,x[iout+j]);
      svmul(0.41,x[iout+j],t2);
      rvec_inc(dip,t2);
    }
    if (bYaw)
      for(m=0; (m<DIM); m++) 
	x[iout+3][m] = x[iout+4][m] = 
	  (1-2*DCONS)*x[iout][m]+DCONS*(x[iout+1][m]+x[iout+2][m]);
  }
  
  box[XX] = 2*cx;
  box[YY] = 2*(cy2+cy);
  box[ZZ] = 2*(1+cz);
  for(i=0; (i<DIM); i++)
    box[i] *= odist;
    
  printf("Unitcell:  %10.5f  %10.5f  %10.5f\n",box[XX],box[YY],box[ZZ]);
  printf("Dipole:    %10.5f  %10.5f  %10.5f (e nm)\n",dip[XX],dip[YY],dip[ZZ]);
}
开发者ID:daniellandau,项目名称:gromacs,代码行数:67,代码来源:mkice.c


示例12: virial

void virial(FILE *fp,gmx_bool bFull,int nmol,rvec x[],matrix box,real rcut,
	    gmx_bool bYaw,real q[],gmx_bool bLJ)
{
  int  i,j,im,jm,natmol,ik,jk,m,ninter;
  rvec dx,f,ftot,dvir,vir,pres,xcmi,xcmj,*force;
  real dx6,dx2,dx1,fscal,c6,c12,vcoul,v12,v6,vctot,v12tot,v6tot;
  t_pbc pbc;
  
  set_pbc(&pbc,box);
  fprintf(fp,"%3s   -  %3s: %6s %6s %6s  %6s %8s %8s %8s\n",
	  "ai","aj","dx","dy","dz","|d|","virx","viry","virz");
  clear_rvec(ftot);
  clear_rvec(vir);
  ninter = 0;
  vctot  = 0;
  v12tot = 0;
  v6tot  = 0;
  natmol = bYaw ? 5 : 3;
  snew(force,nmol*natmol);
  
  for(i=0; (i<nmol); i++) {
    im = natmol*i;
    /* Center of geometry */
    clear_rvec(xcmi);
    for(ik=0; (ik<natmol); ik++)
      rvec_inc(xcmi,x[im+ik]);
    for(m=0; (m<DIM); m++)
      xcmi[m] /= natmol;

    for(j=i+1; (j<nmol); j++) {
      jm = natmol*j;
      /* Center of geometry */
      clear_rvec(xcmj);
      for(jk=0; (jk<natmol); jk++)
	rvec_inc(xcmj,x[jm+jk]);
      for(m=0; (m<DIM); m++)
	xcmj[m] /= natmol;

      /* First check COM-COM distance */
      pbc_dx(&pbc,xcmi,xcmj,dx);
      if (norm(dx) < rcut) {
	ninter++;
	/* Neirest neighbour molecules! */
	clear_rvec(dvir);
	for(ik=0; (ik<natmol); ik++) {
	  for(jk=0; (jk<natmol); jk++) {
	    pbc_dx(&pbc,x[im+ik],x[jm+jk],dx);
	    dx2    = iprod(dx,dx);
	    dx1    = sqrt(dx2);
	    vcoul  = q[ik]*q[jk]*ONE_4PI_EPS0/dx1;
	    vctot += vcoul;
	    
	    if (bLJ) {
	      if (bYaw) {
		c6  = yaw_lj[ik][2*jk];
		c12 = yaw_lj[ik][2*jk+1];
	      }
	      else {
		c6  = spc_lj[ik][2*jk];
		c12 = spc_lj[ik][2*jk+1];
	      }
	      dx6    = dx2*dx2*dx2;
	      v6     = c6/dx6;
	      v12    = c12/(dx6*dx6);
	      v6tot -= v6;
	      v12tot+= v12;
	      fscal  = (vcoul+12*v12-6*v6)/dx2;
	    }
	    else
	      fscal  = vcoul/dx2;
	    for(m=0; (m<DIM); m++) {
	      f[m]     = dx[m]*fscal;
	      dvir[m] -= 0.5*dx[m]*f[m];
	    }
	    rvec_inc(force[ik+im],f);
	    rvec_dec(force[jk+jm],f);
	    /*if (bFull)
	      fprintf(fp,"%3s%4d-%3s%4d: %6.3f %6.3f %6.3f %6.3f"
		      " %8.3f %8.3f %8.3f\n",
		      watname[ik],im+ik,watname[jk],jm+jk,
		      dx[XX],dx[YY],dx[ZZ],norm(dx),
		      dvir[XX],dvir[YY],dvir[ZZ]);*/
	  }
	}
	if (bFull)
	  fprintf(fp,"%3s%4d-%3s%4d: "
		  " %8.3f %8.3f %8.3f\n",
		  "SOL",i,"SOL",j,dvir[XX],dvir[YY],dvir[ZZ]);
	rvec_inc(vir,dvir);
      }
    }
  }
  fprintf(fp,"There were %d interactions between the %d molecules (%.2f %%)\n",
	  ninter,nmol,(real)ninter/(0.5*nmol*(nmol-1)));
  fprintf(fp,"Vcoul: %10.4e  V12: %10.4e  V6: %10.4e  Vtot: %10.4e (kJ/mol)\n",
	  vctot/nmol,v12tot/nmol,v6tot/nmol,(vctot+v12tot+v6tot)/nmol);
  pr_rvec(fp,0,"vir ",vir,DIM,TRUE);
  
  for(m=0; (m<DIM); m++) 
    pres[m] = -2*PRESFAC/(det(box))*vir[m];
//.........这里部分代码省略.........
开发者ID:daniellandau,项目名称:gromacs,代码行数:101,代码来源:mkice.c


示例13: random_h_coords

void random_h_coords(int natmol,int nmol,rvec x[],rvec box,
		     gmx_bool bYaw,real odist,real hdist)
{
#define cx  0.81649658
#define cy  0.47140452
#define cy2 0.94280904
#define cz  0.33333333
  
  rvec xx[24] = {
    { 0,   0,         0 }, /* O1 */
    { 0,   0,         1 }, /* H relative to Oxygen */
    { cx, cy,       -cz },
    { cx, cy,       -cz }, /* O2 */
    { 0, 0,       -1    }, /* H relative to Oxygen */
    { cx,-cy,       +cz },
    { cx, cy+cy2,     0 }, /* O3 */
    { -cx, cy,    -cz   }, /* H relative to Oxygen */
    { 0,   -cy2,    -cz },
    { 0,  2*cy+cy2, -cz }, /* O4 */
    {-cx,-cy,       +cz }, /* H relative to Oxygen */
    { 0 , cy2,      +cz },
    { 0,   0,         1 }, /* O5 */
    {-cx, cy,       +cz }, /* H relative to Oxygen */
    { 0 , -cy2,     +cz },
    { cx, cy,      1+cz }, /* O6 */
    { -cx, -cy,     -cz }, /* H relative to Oxygen */
    { 0,   cy2,     -cz },
    { cx, cy+cy2,     1 }, /* O7 */
    { 0,  0,       -1   }, /* H relative to Oxygen */
    { cx, cy,       +cz },
    { 0,  2*cy+cy2,1+cz }, /* O8 */
    { 0,   0,         1 }, /* H relative to Oxygen */
    { cx,   -cy,    -cz }
  };
  int  i,iin,iout,j,m;
  rvec tmp,t2,dip;
  
  clear_rvec(dip);
  for(i=0; (i<nmol); i++) {
    iin = natmol*i;
    iout = iin;
    svmul(odist,x[iin],x[iout]);
    svmul(-0.82,x[iout],t2);
    rvec_inc(dip,t2);
    for(j=1; (j<=2); j++) {
      svmul(hdist,xx[3*(i % 8)+j],tmp);
      rvec_add(x[iout],tmp,x[iout+j]);
      svmul(0.41,x[iout+j],t2);
      rvec_inc(dip,t2);
    }
  }
  
  box[XX] = 2*cx;
  box[YY] = 2*(cy2+cy);
  box[ZZ] = 2*(1+cz);
  for(i=0; (i<DIM); i++)
    box[i] *= odist;
    
  printf("Unitcell:  %10.5f  %10.5f  %10.5f\n",box[XX],box[YY],box[ZZ]);
  printf("Dipole:    %10.5f  %10.5f  %10.5f (e nm)\n",dip[XX],dip[YY],dip[ZZ]);
}
开发者ID:daniellandau,项目名称:gromacs,代码行数:61,代码来源:mkice.c


示例14: ewald_LRcorrection


//.........这里部分代码省略.........
            if (calc_excl_corr)
            {
                i1  = excl->index[i];
                i2  = excl->index[i+1];

                /* Loop over excluded neighbours */
                for (j = i1; (j < i2); j++)
                {
                    k = AA[j];
                    /*
                     * First we must test whether k <> i, and then, because the
                     * exclusions are all listed twice i->k and k->i we must select
                     * just one of the two.
                     * As a minor optimization we only compute forces when the charges
                     * are non-zero.
                     */
                    if (k > i)
                    {
                        qqA = qiA*chargeA[k];
                        if (qqA != 0.0)
                        {
                            rvec_sub(x[i], x[k], dx);
                            if (bMolPBC)
                            {
                                /* Cheap pbc_dx, assume excluded pairs are at short distance. */
                                for (m = DIM-1; (m >= 0); m--)
                                {
                                    if (dx[m] > 0.5*box[m][m])
                                    {
                                        rvec_dec(dx, box[m]);
                                    }
                                    else if (dx[m] < -0.5*box[m][m])
                                    {
                                        rvec_inc(dx, box[m]);
                                    }
                                }
                            }
                            dr2 = norm2(dx);
                            /* Distance between two excluded particles may be zero in the
                             * case of shells
                             */
                            if (dr2 != 0)
                            {
                                rinv              = gmx_invsqrt(dr2);
                                rinv2             = rinv*rinv;
                                dr                = 1.0/rinv;
#ifdef TABLES
                                r1t               = tabscale*dr;
                                n0                = r1t;
                                assert(n0 >= 3);
                                n1                = 12*n0;
                                eps               = r1t-n0;
                                eps2              = eps*eps;
                                nnn               = n1;
                                Y                 = VFtab[nnn];
                                F                 = VFtab[nnn+1];
                                Geps              = eps*VFtab[nnn+2];
                                Heps2             = eps2*VFtab[nnn+3];
                                Fp                = F+Geps+Heps2;
                                VV                = Y+eps*Fp;
                                FF                = Fp+Geps+2.0*Heps2;
                                vc                = qqA*(rinv-VV);
                                fijC              = qqA*FF;
                                Vexcl            += vc;

                                fscal             = vc*rinv2+fijC*tabscale*rinv;
开发者ID:yhalcyon,项目名称:Gromacs,代码行数:67,代码来源:ewald_util.c


示例15: dielectric


//.........这里部分代码省略.........
        }
        if (time[nfr] <= efit)
        {
            ei = nfr;
        }

        if (bNoJump)
        {

            if (xp)
            {
                remove_jump(fr.box, fr.natoms, xp, fr.x);
            }
            else
            {
                snew(xp, fr.natoms);
            }

            for (i = 0; i < fr.natoms; i++)
            {
                copy_rvec(fr.x[i], xp[i]);
            }

        }

        gmx_rmpbc_trxfr(gpbc, &fr);

        calc_mj(top, ePBC, fr.box, bNoJump, nmols, indexm, fr.x, mtrans[nfr], mass2, qmol);

        for (i = 0; i < isize; i++)
        {
            j = index0[i];
            svmul(top.atoms.atom[j].q, fr.x[j], fr.x[j]);
            rvec_inc(mu[nfr], fr.x[j]);
        }

        /*if(mod(nfr,nshift)==0){*/
        if (nfr%nshift == 0)
        {
            for (j = nfr; j >= 0; j--)
            {
                rvec_sub(mtrans[nfr], mtrans[j], tmp);
                dsp2[nfr-j]  += norm2(tmp);
                xshfr[nfr-j] += 1.0;
            }
        }

        if (fr.bV)
        {
            if (nvfr >= valloc)
            {
                valloc += 100;
                srenew(vfr, valloc);
                if (bINT)
                {
                    srenew(djc, valloc);
                }
                srenew(v0, valloc);
                if (bACF)
                {
                    srenew(cacf, valloc);
                }
            }
            if (time[nfr] <= bvit)
            {
                ii = nvfr;
开发者ID:pjohansson,项目名称:gromacs,代码行数:67,代码来源:gmx_current.cpp


示例16: dd_pmeredist_f

void dd_pmeredist_f(struct gmx_pme_t *pme, pme_atomcomm_t *atc,
                    int n, rvec *f,
                    gmx_bool bAddF)
{
    int *commnode, *buf_index;
    int  nnodes_comm, local_pos, buf_pos, i, scount, rcount, node;

    commnode  = atc->node_dest;
    buf_index = atc->buf_index;

    nnodes_comm = min(2*atc->maxshift, atc->nslab-1);

    local_pos = atc->count[atc->nodeid];
    buf_pos   = 0;
    for (i = 0; i < nnodes_comm; i++)
    {
        scount = atc->rcount[i];
        rcount = atc->count[commnode[i]];
        if (scount > 0 || rcount > 0)
        {
            /* Communicate the forces */
            pme_dd_sendrecv(atc, TRUE, i,
                            atc->f[local_pos], scount*sizeof(rvec),
                            pme->bufv[buf_pos], rcount*sizeof(rvec));
            local_pos += scount;
        }
        buf_index[commnode[i]] = buf_pos;
        buf_pos               += rcount;
    }

    local_pos = 0;
    if (bAddF)
    {
        for (i = 0; i < n; i++)
        {
            node = atc->pd[i];
            if (node == atc->nodeid)
            {
                /* Add from the local force array */
                rvec_inc(f[i], atc->f[local_pos]);
                local_pos++;
            }
            else
            {
                /* Add from the receive buffer */
                rvec_inc(f[i], pme->bufv[buf_index[node]]);
                buf_index[node]++;
            }
        }
    }
    else
    {
        for (i = 0; i < n; i++)
        {
            node = atc->pd[i];
            if (node == atc->nodeid)
            {
                /* Copy from the local force array */
                copy_rvec(atc->f[local_pos], f[i]);
                local_pos++;
            }
            else
            {
                /* Copy from the receive buffer */
                copy_rvec(pme->bufv[buf_index[node]], f[i]);
                buf_index[node]++;
            }
        }
    }
}
开发者ID:FoldingAtHome,项目名称:gromacs,代码行数:70,代码来源:pme-redistribute.c


示例17: calc_order


//.........这里部分代码省略.........
       in index*/
#endif

    teller = 0;

    gpbc = gmx_rmpbc_init(&top->idef, ePBC, natoms);
    /*********** Start processing trajectory ***********/
    do
    {
        if (bSliced)
        {
            *slWidth = box[axis][axis]/nslices;
        }
        teller++;

        set_pbc(&pbc, ePBC, box);
        gmx_rmpbc_copy(gpbc, natoms, box, x0, x1);

        /* Now loop over all groups. There are ngrps groups, the order parameter can
           be calculated for grp 1 to grp ngrps - 1. For each group, loop over all
           atoms in group, which is index[i] to (index[i+1] - 1) See block.h. Of
           course, in this case index[i+1] -index[i] has to be the same for all
           groups, namely the number of tails. i just runs over all atoms in a tail,
           so for DPPC ngrps = 16 and i runs from 1 to 14, including 14
         */


        if (radial)
        {
            /*center-of-mass determination*/
            com[XX] = 0.0; com[YY] = 0.0; com[ZZ] = 0.0;
            for (j = 0; j < comsize; j++)
            {
                rvec_inc(com, x1[comidx[j]]);
            }
            svmul(1.0/comsize, com, com);
        }
        if (distcalc)
        {
            dref[XX] = 0.0; dref[YY] = 0.0; dref[ZZ] = 0.0;
            for (j = 0; j < distsize; j++)
            {
                rvec_inc(dist, x1[distidx[j]]);
            }
            svmul(1.0/distsize, dref, dref);
            if (radial)
            {
                pbc_dx(&pbc, dref, com, dvec);
                unitv(dvec, dvec);
            }
        }

        for (i = 1; i < ngrps - 1; i++)
        {
            clear_rvec(frameorder);

            size = index[i+1] - index[i];
            if (size != nr_tails)
            {
                gmx_fatal(FARGS, "grp %d does not have same number of"
                          " elements as grp 1\n", i);
            }

            for (j = 0; j < size; j++)
            {
                if (radial)
开发者ID:tanigawa,项目名称:gromacs,代码行数:67,代码来源:gmx_order.cpp


示例18: gmx_gyrate


//.........这里部分代码省略.........
        xvgr_legend(out, NLEG, legI, oenv);
    }
    else
    {
        if (bRot)
        {
            if (output_env_get_print_xvgr_codes(oenv))
            {
                fprintf(out, "@ subtitle \"Axes are principal component axes\"\n");
            }
        }
        xvgr_legend(out, NLEG, leg, oenv);
    }
    if (nz == 0)
    {
        gpbc = gmx_rmpbc_init(&top.idef, ePBC, natoms);
    }
    do
    {
        if (nz == 0)
        {
            gmx_rmpbc_copy(gpbc, natoms, box, x, x_s);
        }
        gyro = 0;
        clear_rvec(gvec);
        clear_rvec(gvec1);
        clear_rvec(d);
        clear_rvec(d1);
        for (mol = 0; mol < nmol; mol++)
        {
            tm    = sub_xcm(nz == 0 ? x_s : x, nam, index+mol*nam, top.atoms.atom, xcm, bQ);
            if (nz == 0)
            {
                gyro += calc_gyro(x_s, nam, index+mol*nam, top.atoms.atom,
                                  tm, gvec1, d1, bQ, bRot, bMOI, trans);
            }
            else
            {
                calc_gyro_z(x, box, nam, index+mol*nam, top.atoms.atom, nz, t, out);
            }
            rvec_inc(gvec, gvec1);
            rvec_inc(d, d1);
        }
        if (nmol > 0)
        {
            gyro /= nmol;
            svmul(1.0/nmol, gvec, gvec);
            svmul(1.0/nmol, d, d);
        }

        if (nz == 0)
        {
            if (bRot)
            {
                if (j >= max_moi)
                {
                    max_moi += delta_moi;
                    for (m = 0; (m < DIM); m++)
                    {
                        srenew(moi_trans[m], max_moi*DIM);
                    }
                }
                for (m = 0; (m < DIM); m++)
                {
                    copy_rvec(trans[m], moi_trans[m]+DIM*j);
                }
                fprintf(out, "%10g  %10g  %10g  %10g  %10g\n",
                        t, gyro, d[XX], d[YY], d[ZZ]);
            }
            else
            {
                fprintf(out, "%10g  %10g  %10g  %10g  %10g\n",
                        t, gyro, gvec[XX], gvec[YY], gvec[ZZ]);
            }
        }
        j++;
    }
    while (read_next_x(oenv, status, &t, x, box));
    close_trj(status);
    if (nz == 0)
    {
        gmx_rmpbc_done(gpbc);
    }

    xvgrclose(out);

    if (bACF)
    {
        int mode = eacVector;

        do_autocorr(opt2fn("-acf", NFILE, fnm), oenv,
                    "Moment of inertia vector ACF",
                    j, 3, moi_trans, (t-t0)/j, mode, FALSE);
        do_view(oenv, opt2fn("-acf", NFILE, fnm), "-nxy");
    }

    do_view(oenv, ftp2fn(efXVG, NFILE, fnm), "-nxy");

    return 0;
}
开发者ID:tanigawa,项目名称:gromacs,代码行数:101,代码来源:gmx_gyrate.cpp


示例19: calc_chargegroup_radii

void calc_chargegroup_radii(const gmx_mtop_t *mtop,rvec *x,
                            real *rvdw1,real *rvdw2,
                            real *rcoul1,real *rcoul2)
{
    real r2v1,r2v2,r2c1,r2c2,r2;
    int  ntype,i,j,mb,m,cg,a_mol,a0,a1,a;
    gmx_bool *bLJ;
    gmx_molblock_t *molb;
    gmx_moltype_t *molt;
    t_block *cgs;
    t_atom *atom;
    rvec cen;

    r2v1 = 0;
    r2v2 = 0;
    r2c1 = 0;
    r2c2 = 0;

    ntype = mtop->ffparams.atnr;
    snew(bLJ,ntype);
    for(i=0; i<ntype; i++)
    {
        bLJ[i] = FALSE;
        for(j=0; j<ntype; j++)
        {
            if (mtop->ffparams.iparams[i*ntype+j].lj.c6  != 0 ||
                mtop->ffparams.iparams[i*ntype+j].lj.c12 != 0)
            {
                bLJ[i] = TRUE;
            }
        }
    }

    a_mol = 0;
    for(mb=0; mb<mtop->nmolblock; mb++)
    {
        molb = &mtop->molblock[mb];
        molt = &mtop->moltype[molb->type];
        cgs  = &molt->cgs;
        atom = molt->atoms.atom;
        for(m=0; m<molb->nmol; m++)
        {
            for(cg=0;  

鲜花

握手

雷人

路过

鸡蛋
该文章已有0人参与评论

请发表评论

全部评论

专题导读
上一篇:
C++ rvmExceptionCheck函数代码示例发布时间:2022-05-30
下一篇:
C++ rvec_dec函数代码示例发布时间:2022-05-30
热门推荐
阅读排行榜

扫描微信二维码

查看手机版网站

随时了解更新最新资讯

139-2527-9053

在线客服(服务时间 9:00~18:00)

在线QQ客服
地址:深圳市南山区西丽大学城创智工业园
电邮:jeky_zhao#qq.com
移动电话:139-2527-9053

Powered by 互联科技 X3.4© 2001-2213 极客世界.|Sitemap