LAMMPS WWW Site - LAMMPS Documentation - LAMMPS Mailing List Archives
Re: [lammps-users] LAMMPS library: Gather specific atoms and insert small molecules
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: [lammps-users] LAMMPS library: Gather specific atoms and insert small molecules


From: Axel Kohlmeyer <akohlmey@...24...>
Date: Sat, 2 Sep 2017 10:17:59 -0400



On Fri, Sep 1, 2017 at 7:44 PM, Stephan Grein <stephan.grein@...5733...> wrote:
Dear Axel,

thanks, this sounds good.

What I'm not understanding so far is,
that I thought by using something like

double *coords = (double*)lammps_extract_variable(lmp,"MYCOORDS","MYGRP");

I would get only the coordinates for the atoms in MYGRP
and thus would not have to use gmask as selector?

​no. atom style variables are always returned as​ an array of the size of the number of local atoms.
 


A last question (hopefully) concerning this:
If I run that code now with MPI (more than one proc),
I get some garbage output, like:
x[1005]=7199496020021154586587594...

Is this to be expected without using the selector?

​no. this is because the demo code i wrote only works for 1 processor.  for 1 processor nlocal == natoms.
atom style variables are stored distributed across the domains, so the length of the array is nlocal and not natoms for more than one MPI rank.

sorry about that. so the more correct library interface test code looks like this:

#include <mpi.h>
#include <stdio.h>
#include <inttypes.h>
#include "library.h"

int main (int argc, char **argv)

{
    void *lmp;
    int i,gatoms;
    int nlocal;

    MPI_Init(&argc,&argv);
    lammps_open(0,NULL,MPI_COMM_WORLD,&lmp);
    lammps_file(lmp,"in.protein");
    lammps_command(lmp,"group MYGRP type 2");
    lammps_command(lmp,"variable MYCOORDS atom x");
    lammps_command(lmp,"variable MYIDS atom id");
    lammps_command(lmp,"variable MYCOUNT equal count(MYGRP)");

    nlocal = *(int *)lammps_extract_global(lmp,"nlocal");
    gatoms = (int) *(double *)lammps_extract_variable(lmp,"MYCOUNT",NULL);
    printf("local atoms = %d  atoms in group = %d\n",nlocal,gatoms);
    double *coords = (double *)lammps_extract_variable(lmp,"MYCOORDS","MYGRP");
    double *ids = (double *)lammps_extract_variable(lmp,"MYIDS","MYGRP");
    for (i=0; i < nlocal; ++i) {
        if (ids[i] != 0.0)
            printf("x[%03d]=%10.4f\n",(int)ids[i],coords[i]);
    }
    lammps_free((void *)coords);
    lammps_close(lmp);
    MPI_Finalize();
    return 0;
}



​axel.​

 

Thanks in advance.

Best regards,
Stephan

Am 02/09/17 um 01:36 schrieb Axel Kohlmeyer:
>
>
> On Fri, Sep 1, 2017 at 7:14 PM, Stephan Grein
> <stephan.grein@...7108.....frankfurt.de
> <mailto:stephan.grein@...7107...uni-frankfurt.de>> wrote:
>
>     Dear Axel,
>
>     thanks for your excellent help!
>
>     I used gatoms in the for loop, not natoms.
>     Now the question is, how could I get only
>     the x coord of the three atoms in the group
>     without checking for inequality with 0.0?
>
>
> ​you could pull in another atom style variable that contains the output
> of ​gmask(GROUP), which is 1 for group members and 0 for other atoms,
> and then use this variable as selector.
>
> axel.
>
>  
>
>
>     Best regards,
>     Stephan
>
>
>
> --
> Dr. Axel Kohlmeyer  akohlmey@...24... <mailto:akohlmey@...24...>
>  http://goo.gl/1wk0
> College of Science & Technology, Temple University, Philadelphia PA, USA
> International Centre for Theoretical Physics, Trieste. Italy.




--
Dr. Axel Kohlmeyer  akohlmey@...24...  http://goo.gl/1wk0
College of Science & Technology, Temple University, Philadelphia PA, USA
International Centre for Theoretical Physics, Trieste. Italy.