Problems with MEAM library for couple code

Hello,

My group and I are using LAMMPS as a library for more than 2 years now. LAMMPS is coupled to our Fortran code (kinetic Activation-Relaxation Technique kART) just to compute forces and use the wide catalog of implemented potentials.

We found something that looks like a bug. let me explain :

we have two instance pointers one for global force calculations and a second for local force calculations.
So that we have two lammps_open() command with pointer named lmp and lmp_local
and also two times the lammps_file() command to read the same input file but to two different instances.

Everything works well with build in potentials as SW and LJ it is also working with reaxff potential.
Our code decide if we want to compute force locally by calling lmp_local or to compute force globally by calling lmp both instance are open on the same core so that the two lammps instances are open and waiting for instructions from our code kART.

Now, with the MEAM potential, we are not able to pass throughout the initialization of LAMMPS because during the second lammps_file() the code stops. We have done some work and it seems that reading a second time the pair_coeff line is the source of the problem. But the thing is that if we remove this line the lmp_local instance are not able to compute force because we did not provide the potential parameter.

The only difference in this line compare to other potential is that you can pass two files to the pair_coeff instead of one usually.

See here the error message that we get with our code

forrtl: severe (151): allocatable array is already allocated
Image PC Routine Line Source
libintlc.so.5 00002AB8E4C71F7A Unknown Unknown Unknown
libintlc.so.5 00002AB8E4C70AF5 Unknown Unknown Unknown
libifcore.so.5 00002AB8E44311F2 Unknown Unknown Unknown
libifcore.so.5 00002AB8E43A65FB Unknown Unknown Unknown
libifcore.so.5 00002AB8E43F37A3 Unknown Unknown Unknown
KMCART_briareelam 0000000000DD6CF6 compute_pair_meam 187 meam_setup_done.F
KMCART_briareelam 0000000000DD5E76 meam_setup_done_ 84 meam_setup_done.F
KMCART_briareelam 0000000000C78DC6 _ZN9LAMMPS_NS8Pai 369 pair_meam.cpp
KMCART_briareelam 0000000000B0506E ZN9LAMMPS_NS5Inp 1626 input.cpp
KMCART_briareelam 0000000000AFFACD ZN9LAMMPS_NS5Inp 675 input.cpp
KMCART_briareelam 0000000000AFD9DF ZN9LAMMPS_NS5Inp 243 input.cpp
KMCART_briareelam 0000000000AFDC22 ZN9LAMMPS_NS5Inp 277 input.cpp
KMCART_briareelam 00000000008E3A11 lammps_file 100 library.cpp
KMCART_briareelam 000000000077746E lammps_mp_lammps
302 LAMMPS.f90
KMCART_briareelam 00000000005B2030 init_potential_la 119 calcfo_lammps.f90
KMCART_briareelam 00000000005BFD9E calcforce_mod_mp
175 calcforce.f90
KMCART_briareelam 00000000007BFFBE MAIN
139 main_KART.f90
KMCART_briareelam 00000000004A0B9C Unknown Unknown Unknown
libc.so.6 000000332461ED5D Unknown Unknown Unknown
KMCART_briareelam 00000000004A0A99 Unknown Unknown Unknown

Hello,

My group and I are using LAMMPS as a library for more than 2 years now.
LAMMPS is coupled to our Fortran code (kinetic Activation-Relaxation
Technique kART) just to compute forces and use the wide catalog of
implemented potentials.

We found something that looks like a bug. let me explain :

we have two instance pointers one for global force calculations and a second
for local force calculations.
So that we have two lammps_open() command with pointer named lmp and
lmp_local
and also two times the lammps_file() command to read the same input file but
to two different instances.

i don't think this can work at all with MEAM, since that is
implemented as fortran library with no support for multiple instances.

axel.

All the C++/C code in LAMMPS, including all the external libs

under the lib directory (so far as I know) is written in a manner

that there are no global dependencies. So you can instantiate

multiple instances of LAMMPS. However as Axel says, the

MEAM lib was written in Fortran and it may well be that

there are effectively global variables which would mean

it cannot be instantiated twice. Whether that is always the

case with a Fortran lib (doubt it), or just the way that the

MEAM lib was written, I’m not sure. If the latter, there is

probably a way to change that if someone better at modern

Fortran knows what causes those kind of issues.

Steve

Hi Steve and Axel,

Thank you for your answers.

I manage to solve the problem by modifying the MEAM library.
The fix is simple but you might not like it.
I just put if statements as follow in team_setup_done.F

c allocate memory for array that defines the potential
if (not(allocated(phir))) then
allocate(phir(nr,(neltypes*(neltypes+1))/2))
endif
c allocate coeff memory

if (not(allocated(phirar))) then
allocate(phirar(nr,(neltypes*(neltypes+1))/2))
endif
if (not(allocated(phirar1))) then
allocate(phirar1(nr,(neltypes*(neltypes+1))/2))
endif
if (not(allocated(phirar2))) then
allocate(phirar2(nr,(neltypes*(neltypes+1))/2))
endif
if (not(allocated(phirar3))) then
allocate(phirar3(nr,(neltypes*(neltypes+1))/2))
endif
if (not(allocated(phirar4))) then
allocate(phirar4(nr,(neltypes*(neltypes+1))/2))
endif
if (not(allocated(phirar5))) then
allocate(phirar5(nr,(neltypes*(neltypes+1))/2))
endif
if (not(allocated(phirar6))) then
allocate(phirar6(nr,(neltypes*(neltypes+1))/2))
endif

So that if these arrays are already allocate we are not trying to reallocated them.

This is good for us since we want the same set of parameter read twice, but it might not be a good fix for everyone because one might want to read two different sets of MEAM parameters.

PS : I think that our two-three years of work on the kART (a kinetic Monte Carlo method) LAMMPS coupling (with LAMMPS as a library) is good enough to ask you, if we can be add to your lammps Project in the COUPLE folder.
Please tell me, what are the requirement to do so.

See here two references of our work (all result in my article calculate force using Stillinger-Weber potential using LAMMPS) :

Mickaël Trochet, et al##### Phys. Rev. B 91, 224106 – Published 16 June 2015Normand Mousseau, et al

Comput. Mater. Sci. 100,111 (2015)

Mickaël Trochet, PhD student

Département de Physique A-408

Université de Montréal
2900 blvd. Édouard-Monpetit, Montréal, QC, Canada, H3T 1J4

email: mickaeltrochet@…472…, mickael.laurent.trochet@…1077…

All the C++/C code in LAMMPS, including all the external libs

under the lib directory (so far as I know) is written in a manner

that there are no global dependencies. So you can instantiate

multiple instances of LAMMPS. However as Axel says, the

MEAM lib was written in Fortran and it may well be that

there are effectively global variables which would mean

it cannot be instantiated twice. Whether that is always the

case with a Fortran lib (doubt it), or just the way that the

MEAM lib was written, I’m not sure. If the latter, there is

probably a way to change that if someone better at modern

Fortran knows what causes those kind of issues.

Steve

Hi Steve and Axel,

Thank you for your answers.

I manage to solve the problem by modifying the MEAM library.
The fix is simple but you might not like it.
I just put if statements as follow in team_setup_done.F

at this point it might be worth having another look at the extension
to the meam/spline code that i have sitting in my inbox for almost two
years now without ever finding the time to finish the integration into
LAMMPS.

these kind of things usually require somebody that has a vested
interest to test the code after it has been merged.
the upshot would be, that the resulting simulations might be
significantly faster (at least as far as LAMMPS is concerned) and the
limitations of the fortran library would be less relevant as well.

axel.

You’re correct that what you’ve done is not a general solution b/c

it doesn’t allow the lib to use different params. But it’s good

you identified the problem. Those arrays, along with many

other variables are defined in mean_data.F. Does the

way they are defined mean they are effectively global variables,

owned in common by many instances of the lib?

If so, is there some way in Fortran to make that module data

be private to an instance and not global? So that each instance

has its own copy of the data? That would likely solve all the

problems.

Re: adding kART to the COUPLE folder, that sounds like a good

complex example to have, assuming it is not a large volume

of code/files. See this section of the manual for general guidelines

on submitting content for the code or other tools:

http://lammps.sandia.gov/doc/Section_modify.html#submitting-new-features-for-inclusion-in-lammps

Steve

All FORTRAN variables are by default PUBLIC (ie global) unless declared PRIVATE to the declaring MODULE. The same applies to FUNCTIONS which are declared within a MODULE.

You're correct that what you've done is not a general solution b/c
it doesn't allow the lib to use different params. But it's good
you identified the problem. Those arrays, along with many
other variables are defined in mean_data.F. Does the
way they are defined mean they are effectively global variables,
owned in common by many instances of the lib?

yes.

If so, is there some way in Fortran to make that module data
be private to an instance and not global? So that each instance
has its own copy of the data? That would likely solve all the
problems.

in principle, yes, but it would be complicated.
the most straightforward way to deal with this would be use storage
allocated by the pair style (i.e. in C++) for this.
otherwise you'd have to create a derived type (like a struct in
C/C++) and pass a reference to an instance of the derived type back up
the the pair style and manage it from there and then deallocate it
from the destructor. basically, you'd have to treat it similar to how
the C library interface handles LAMMPS instances.

axel.

I just realized that the reaxff fortran library allows multiple instance and read of files.

I think that it can be a good starting point to see how it is done in reaxff to have an idea of a fix for meam fortran library.

Also our code is written in fortran and (base on my short experience 4 years of PhD ) I think that it absolutely possible to create local variable that will set and allocate arrays of data related to potential parameter, for each instance separately. The problem is that I don’t know how much time and how many modifications to the meam files it will take.

For adding kART to the COUPLE example, I spoke with my supervisor about it and he thinks that we need to work more on the setup simplicity and also a good clean up before bothering you with it.
But I hope that I will be able to do that before the end of my PhD, because I think it will be a good achievement about my computation development part and a good addition to both code.

Sincerely Mickaël,

Mickaël Trochet, PhD student

Département de Physique A-408

Université de Montréal
2900 blvd. Édouard-Monpetit, Montréal, QC, Canada, H3T 1J4

email: mickaeltrochet@…472…, mickael.laurent.trochet@…1077…

You’re correct that what you’ve done is not a general solution b/c

it doesn’t allow the lib to use different params. But it’s good

you identified the problem. Those arrays, along with many

other variables are defined in mean_data.F. Does the

way they are defined mean they are effectively global variables,

owned in common by many instances of the lib?

If so, is there some way in Fortran to make that module data

be private to an instance and not global? So that each instance

has its own copy of the data? That would likely solve all the

problems.

Re: adding kART to the COUPLE folder, that sounds like a good

complex example to have, assuming it is not a large volume

of code/files. See this section of the manual for general guidelines

on submitting content for the code or other tools:

http://lammps.sandia.gov/doc/Section_modify.html#submitting-new-features-for-inclusion-in-lammps

Steve

All FORTRAN variables are by default PUBLIC (ie global) unless declared
PRIVATE to the declaring MODULE. The same applies to FUNCTIONS which are
declared within a MODULE.

public vs. private is only managing read/write access. (same as
public/private in c++ classes), but you still can have only one
instance of a module.

in order to have multiple instances of an entity you need a way to
associate a specific instance of the storage locations with the
corresponding instance of the C++ class.
there are several ways of doing it. i just mentioned two. i can think
of at least two more (that are more complex or ugly).

axel.

I just realized that the reaxff fortran library allows multiple instance and
read of files.

I think that it can be a good starting point to see how it is done in reaxff
to have an idea of a fix for meam fortran library.

Also our code is written in fortran and (base on my short experience 4 years
of PhD ) I think that it absolutely possible to create local variable that
will set and allocate arrays of data related to potential parameter, for
each instance separately. The problem is that I don't know how much time and
how many modifications to the meam files it will take.

it will be more than you currently think. while you don't get an error
with your hack, there are more places where there is storage that is
dependent on specific input. you just get away with it, because you
use the same parameters for all instances. creating "local" variables
is not sufficient. you have to be able to associate them with a
specific instance of the meam pair style inside each LAMMPS instance.

axel.

I agree with Axel that one solution is let LAMMPS

handle the MEAM memory in the pair_meam.cpp wrapper

of the lib. But I much prefer not to do that, b/c

there’s a lot of data (arrays and scalars) in meam_data.F.

If the MEAM lib were writtten in C/C++, this would

be simple. The lib would allocate everything local each

time it was instantiated.

There must be a way to do something similar in

Fortran with how/when it allocates memory ??

So that the caller (LAMMPS in this case) doesn’t

need to know anything about it. Otherwise there

would be no examples of Fortran libs which can

be used multiple times by a single caller ?

How does the Fortran ReaxFF lib handle its internal

memory? I’m actually surprised that it can in fact

be multiply instantiated by LAMMPS. I didn’t think

that was the case.

Steve

I agree with Axel that one solution is let LAMMPS
handle the MEAM memory in the pair_meam.cpp wrapper
of the lib. But I much prefer not to do that, b/c
there’s a lot of data (arrays and scalars) in meam_data.F.

If the MEAM lib were writtten in C/C++, this would
be simple. The lib would allocate everything local each
time it was instantiated.

A module in fortran is very similar to a namespace in c++.

Using local allocation as such is not enough.

As you need to store those in a location that is unique to that instance of the pair style and can persist between multiple calls to the library.

Axel.

I agree with Axel that one solution is let LAMMPS
handle the MEAM memory in the pair_meam.cpp wrapper
of the lib. But I much prefer not to do that, b/c
there’s a lot of data (arrays and scalars) in meam_data.F.

If the MEAM lib were writtten in C/C++, this would
be simple. The lib would allocate everything local each
time it was instantiated.

There must be a way to do something similar in
Fortran with how/when it allocates memory ??
So that the caller (LAMMPS in this case) doesn’t
need to know anything about it. Otherwise there
would be no examples of Fortran libs which can
be used multiple times by a single caller ?

Almost all of those require the “workspace” storage to be allocated by the caller. C is not different unless the library returns a handle, e.g. the FILE handle of the stdio library.

Axel

Ok so if I understand well, a solution will be to modify the pair_meam.cpp wrapper of LAMMPS. If this is the case, I can’t help you because my c++ skill are close to nothing.

Just to be clear one more time, I don’t know how reaxff handle the potential parameter separately between instance of LAMMPS, but I see only two solutions it is either inside the Library and then written in fortran or in pair_reax.cpp and then written in c++. However, I assure you that it works fine with reaxff.

I’m sorry to not be able to help you more than that.

Sincerely Mickaël,

Mickaël Trochet, PhD student

Département de Physique A-408

Université de Montréal
2900 blvd. Édouard-Monpetit, Montréal, QC, Canada, H3T 1J4

email: mickaeltrochet@…472…, mickael.laurent.trochet@…5968…

I agree with Axel that one solution is let LAMMPS
handle the MEAM memory in the pair_meam.cpp wrapper
of the lib. But I much prefer not to do that, b/c
there’s a lot of data (arrays and scalars) in meam_data.F.

If the MEAM lib were writtten in C/C++, this would
be simple. The lib would allocate everything local each
time it was instantiated.

There must be a way to do something similar in
Fortran with how/when it allocates memory ??
So that the caller (LAMMPS in this case) doesn’t
need to know anything about it. Otherwise there
would be no examples of Fortran libs which can
be used multiple times by a single caller ?
Almost all of those require the “workspace” storage to be allocated by the caller. C is not different unless the library returns a handle, e.g. the FILE handle of the stdio library.

Axel

Ok so if I understand well, a solution will be to modify the pair_meam.cpp
wrapper of LAMMPS. If this is the case, I can't help you because my c++
skill are close to nothing.

Just to be clear one more time, I don't know how reaxff handle the potential
parameter separately between instance of LAMMPS, but I see only two
solutions it is either inside the Library and then written in fortran or in
pair_reax.cpp and then written in c++. However, I assure you that it works
fine with reaxff.

you may be able to run, and you may get lucky because your multiple
instances are effectively describing the same system, if i understood
correctly. this must be the reason that you get away with the crazy
hack you are doing for MEAM, you still *will* be constantly
overwriting storage used by the multiple instances. if i read the REAX
package code correctly, it uses common blocks on the fortran side
(i.e. its programming style predates that of the MEAM library) and
makes them visible to C++ by declaring them global struct variables
with C-style bindings in the C++ code. that is in no way reentrant. in
all honesty, i would be very concerned about this, if i were you. if
it has worked and produced 100% correct results, then only be extreme
chance and you should consider yourself one extremely lucky
individual. i would encourage you to compare your results with what
you would get with reax/c instead of reax.

as it so happens i am this month teaching a couple classes in the HPC
master program in Trieste and i am also touching the very subject of
this matter during this time (not LAMMPS, but Fortran and C++
interoperability and resulting problems). in preparation i have spent
some time looking into this (and worse things), so it is not that i am
making this up as i go along... :wink:

axel.

Hi,

While i have not studied the MEAM or LAMMPS code, my free advice is to
follow the rewrite MEAM in C++ path.
I agree with the comments of Axel in this thread.
My latest experience of several relevant ones was assisting with
debugging the result of converting an Amber Fortran program into a library.
The code was much larger, at least 175000 lines, and legacy issues were
manifold, but it was done by a seasoned programmer.

Based on this size of MEAM analysis, a conversion project might be a
good exercise in learning some LAMMPS development for a graduate student.
Of course not being even a LAMMPS user i have no perspective on the
importance of this.

$ /usr/local/lammps-7Dec15/lib/meam wc *F
   26 36 590 meam_cleanup.F
   86 429 3453 meam_data.F
  282 848 9196 meam_dens_final.F
  564 1801 18306 meam_dens_init.F
  608 1982 23692 meam_force.F
1019 2932 30791 meam_setup_done.F
  111 322 3041 meam_setup_global.F
  204 560 5404 meam_setup_param.F
2900 8910 94473 total

good luck,
scott

Sure, we’d add that to LAMMPS as a meam/c contribution.

I’ll note that the person who did this for ReaxFF, producing reax/c
for LAMMPS along the way, received a PhD for his efforts - he did a lot of other good stuff
as well :slight_smile:

Steve