On Wednesday 30 June 2010 19:45:31 Mikheil Azatov wrote:
> Thank you so much Dr. Axel Arnold for your response. I've been trying to
> make necessary changes to the code but I encountered a lot of questions
> while understanding how all the data parameters work. I would really
> appreciate if you could help me out with some of them.
>
> This is the code that I was thinking to add to the force.c with some
> questions within the code and after.
>
>
> for (c = 0; c < local_cells.n; c++) {
>
> for (n = 0; n < dd.cell_inter[c].n_neighbors; n++) {
> pairs = dd.cell_inter[c].nList[n].vList.pair;
> np = dd.cell_inter[c].nList[n].vList.n;
>
> /* verlet list loop */
> for(i=0; i<2*np; i+=2) {
> p1 = pairs[i]; /* pointer to particle 1 */
> p2 = pairs[i+1]; /* pointer to particle 2 */
> dist2 = distance2vec(p1->r.p, p2->r.p, vec21);
>
> angle = /* angle between paticle 1 particle 2 and third particle
> with id n_id */
> /*Question: I'm not sure how to find the particle with let's say
> identity n_id ? */
> /*Is it just local_particles[n_id]*/
Yes, exactly. Provided, the particle is stored at all on this node, which is
only the case if it is not too far from the other ones. If the particles does
not exist on the local node, local_particles[n_id] will be NULL. In this case,
you anyways cannot specify an angle potential, since also for this potential,
the distance would be too big. The best solution would be to raise a
background error and do nothing. Look in forces.h for runtime_error, how to do
that - that is most similar, since also there, bonds raise errors if the
bonding partners end up on different nodes.
>
> /*now I'm adding the bond if the distance and angle are small
> enough*/
> if (dist2<"bind_distance")&&( angle < "bind_angle" ){
> int part= p1->p.identity; /*id of particle 1*/
> /*Is the value of int part in local_change_bond the same as ID of
> the particle ?*/
Yes, which you can also see from the fact that the local_change_bond uses
local_particles[part] :-).
> int *bond;
That won't work, you need to allocate space here. So, should be "int bond[3]".
> *bond[0]=bond_id(again from espresso file);
> *bond[1]=p1->p.identity;
> *bond[1]=p2->p.identity;
> /*So as I understood bond[0] should have id of the interaction, and
> the next entries should be the particles that are acatually
> interacting with each other*/
Yes, the id of the bond is the number that you gave to the "inter <bond_id>
angle..." command. Best would be to make that a global variable (see below),
so that you can change it in the script.
> int delete=0
> local_change_bond(int part, int *bond, int delete)
> /*Then same for angle bond*/
> /*Question: For angle bond as I understand we will have *bond[1] as
> middle particle *bond[2] left, bond[3], right ?*/
Yes, exactly.
>
> }
> }
> }
>
>
> Questions:
>
> 1)How to access parameters from espresso file? For example let's say I have
> bind_distance(or bind_angle) in Espresso file. How I can I access them here
> ?
That would be done via global variables, which is handled by
global.c/global.h. For example, to get the bond_id as above, you need the
following at global level:
in forces.c:
int bond_id;
/* function to check validity of the provided bond_id */
int bond_id_callback(Tcl_Interp *interp, void *_data)
{
bond_id = *(int *)_data;
if (bond_id < 0 || bond_id >= n_bonded_ia) {
Tcl_AppendResult(interp, "bond_id out of range", (char *)NULL);
return (TCL_ERROR);
}
mpi_bcast_parameter(FIELD_BONDID);
return (TCL_OK);
}
in forces.h:
/* export the static variable in forces.c to the rest of the code */
extern int bond_id;
and in global.c/.h add a new field at the end. That is, in global.h add a
define for the FIELD code, and in global.c, add the descriptor line, which in
this case would be
{&bond_id, TYPE_INT, 1, "bond_id", bond_id_callback, 2 },
The same you do also for the bond_distance etc parameters.
> 2)Should I add this part of code right before the init_forces in force calc
> ?
Doesn't matter at all, you can do both, I think.
> 3)To make new code run, should I compile Espresso with new code in the
> different directory ?
That I don't understand. You have to add the things to the Espresso code, and
then compile Espresso again. Of course, you can make a copy of the source
before, but in any case, you need to modify directly the sources.
> 4)Also it seems to me that this way I compare only the particles from
> neighbor cells of cell c and not within cell c itself. Is it right ? If so
> , how should I check the particles withing the cell ?
No, as far as I remember, the neighboring particles from the same cell are
also in the Verlet list; that is the cell is per definition a neighbor of
itself. However, most probably most of your bonds will be within one cell;
therefore, if that would be different, you would easily see it.
So far, looks like you are progressing fast :-). If there are more questions,
don't hesitate to ask, I try to answer as quick as possible. Just this week
was a bit busy....
Many regards,
Axel
--
JP Dr. Axel Arnold Tel: +49 711 685 67609
ICP, Universität Stuttgart Email:
address@hidden
Pfaffenwaldring 27
70569 Stuttgart, Germany