git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@8132 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp 2012-05-22 15:27:35 +00:00
parent ab865c0b2d
commit 10507a091e
2 changed files with 167 additions and 12 deletions

View File

@ -1,15 +1,21 @@
This directory has a simple C and C++ code that shows how LAMMPS can
be linked to a driver application as a library. The purpose is to
illustrate how another code could perform computations while using
LAMMPS to perform MD, or how an umbrella code or script could call
both LAMMPS and some other code to perform a coupled calculation.
This directory has a simple C, C++, and Fortran code that shows how
LAMMPS can be linked to a driver application as a library. The purpose
is to illustrate how another code could perform computations while
using LAMMPS to perform MD on all or a subset of the processors, or
how an umbrella code or script could call both LAMMPS and some other
code to perform a coupled calculation.
simple.cpp is the C++ driver
simple.c is the C driver
simple.cpp is the C++ driver
simple.c is the C driver
simple.f90 is the Fortran driver
libfwrapper.c is the Fortran-to-C wrapper
The 2 codes do the same thing, so you can compare them to see how to
drive LAMMPS in this manner. The C driver is similar in spirit to
what one could use from a Fortran program or scripting language.
The 3 codes do the same thing, so you can compare them to see how to
drive LAMMPS in this manner. The C driver is similar in spirit to what
one could use to write a scripting language interface. The Fortran
driver in addition requires a wrapper library that interfaces the C
interface of the LAMMPS library to Fortran and also translates the MPI
communicator from Fortran to C.
You can then build either driver code with a compile line something
like this, which includes paths to the LAMMPS library interface, MPI,
@ -28,8 +34,17 @@ gcc -I/home/sjplimp/lammps/src -c simple.c
gcc -L/home/sjplimp/lammps/src simple.o \
-llmp_g++ -lfftw -lmpich -lpthread -lstdc++ -o simpleC
You then run simpleCC or simpleC on a parallel machine on some number
of processors Q with 2 arguments:
This builds the Fortran wrapper and driver with the LAMMPS library
using a Fortran and C compiler:
cp ../fortran/libfwrapper.c .
gcc -I/home/sjplimp/lammps/src -c libfwrapper.c
gfortran -I/home/sjplimp/lammps/src -c simple.f90
gfortran -L/home/sjplimp/lammps/src simple.o libfwrapper.o \
-llmp_g++ -lfftw -lfmpich -lmpich -lpthread -lstdc++ -o simpleF
You then run simpleCC, simpleC, or simpleF on a parallel machine
on some number of processors Q with 2 arguments:
mpirun -np Q simpleCC P in.lj
@ -49,6 +64,11 @@ The C driver is calling C-style routines in the src/library.cpp file
of LAMMPS. You could add any functions you wish to this file to
manipulate LAMMPS data however you wish.
The Fortran driver is using the same C-style routines, but requires an
additional wrapper to make them Fortran callable. Only a subset of the
library functions are currently wrapped, but it should be clear how to
extend the wrapper if desired.
The C++ driver does the same thing, except that it instantiates LAMMPS
as an object first. Some of the functions in src/library.cpp can be
invoked directly as methods within appropriate LAMMPS classes, which

135
couple/simple/simple.f90 Normal file
View File

@ -0,0 +1,135 @@
! LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
! www.cs.sandia.gov/~sjplimp/lammps.html
! Steve Plimpton, sjplimp@sandia.gov, Sandia National Laboratories
!
! Copyright (2003) Sandia Corporation. Under the terms of Contract
! DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
! certain rights in this software. This software is distributed under
! the GNU General Public License.
!
! See the README file in the top-level LAMMPS directory.
! f_driver = simple example of how an umbrella program
! can invoke LAMMPS as a library on some subset of procs
! Syntax: f_driver P in.lammps
! P = # of procs to run LAMMPS on
! must be <= # of procs the driver code itself runs on
! in.lammps = LAMMPS input script
! See README for compilation instructions
PROGRAM f_driver
IMPLICIT NONE
INCLUDE 'mpif.h'
INTEGER, PARAMETER :: fp=20
INTEGER :: n, narg, ierr, me, nprocs, natoms
INTEGER :: lammps, nprocs_lammps, comm_lammps
INTEGER (kind=8) :: ptr
REAL (kind=8), ALLOCATABLE :: x(:)
REAL (kind=8), PARAMETER :: epsilon=0.1
CHARACTER (len=64) :: arg
CHARACTER (len=1024) :: line
! setup MPI and various communicators
! driver runs on all procs in MPI_COMM_WORLD
! comm_lammps only has 1st P procs (could be all or any subset)
CALL mpi_init(ierr)
narg = command_argument_count()
IF (narg /= 2) THEN
PRINT *, 'Syntax: f_driver P in.lammps'
CALL mpi_abort(MPI_COMM_WORLD,1,ierr)
END IF
CALL mpi_comm_rank(MPI_COMM_WORLD,me,ierr);
CALL mpi_comm_size(MPI_COMM_WORLD,nprocs,ierr);
CALL get_command_argument(1,arg)
READ (arg,'(I10)') nprocs_lammps
IF (nprocs_lammps > nprocs) THEN
IF (me == 0) THEN
PRINT *, 'ERROR: LAMMPS cannot use more procs than available'
CALL mpi_abort(MPI_COMM_WORLD,2,ierr)
END IF
END IF
lammps = 0
IF (me < nprocs_lammps) THEN
lammps = 1
ELSE
lammps = MPI_UNDEFINED
END IF
CALL mpi_comm_split(MPI_COMM_WORLD,lammps,0,comm_lammps,ierr)
! open LAMMPS input script on rank zero
CALL get_command_argument(2,arg)
OPEN(UNIT=fp, FILE=arg, ACTION='READ', STATUS='OLD', IOSTAT=ierr)
IF (ierr /= 0) THEN
PRINT *, 'ERROR: Could not open LAMMPS input script'
CALL mpi_abort(MPI_COMM_WORLD,3,ierr);
END IF
! run the input script thru LAMMPS one line at a time until end-of-file
! driver proc 0 reads a line, Bcasts it to all procs
! (could just send it to proc 0 of comm_lammps and let it Bcast)
! all LAMMPS procs call lammps_command() on the line */
IF (lammps == 1) CALL lammps_open(comm_lammps,ptr)
n = 0
DO
IF (me == 0) THEN
READ (UNIT=fp, FMT='(A)', IOSTAT=ierr) line
n = 0
IF (ierr == 0) THEN
n = LEN(TRIM(line))
IF (n == 0 ) THEN
line = ' '
n = 1
END IF
END IF
END IF
CALL mpi_bcast(n,1,MPI_INTEGER,0,MPI_COMM_WORLD,ierr)
IF (n == 0) EXIT
CALL mpi_bcast(line,n,MPI_CHARACTER,0,MPI_COMM_WORLD,ierr)
IF (lammps == 1) CALL lammps_command(ptr,line,n)
END DO
CLOSE(UNIT=fp)
! run 10 more steps
! get coords from LAMMPS
! change coords of 1st atom
! put coords back into LAMMPS
! run a single step with changed coords */
IF (lammps == 1) THEN
CALL lammps_command(ptr,'run 10',6)
CALL lammps_get_natoms(ptr,natoms)
ALLOCATE(x(3*natoms))
CALL lammps_get_coords(ptr,x)
x(1) = x(1) + epsilon
CALL lammps_put_coords(ptr,x)
DEALLOCATE(x)
CALL lammps_command(ptr,'run 1',5);
END IF
! free LAMMPS object
IF (lammps == 1) CALL lammps_close(ptr);
! close down MPI
CALL mpi_finalize(ierr)
END PROGRAM f_driver