diff --git a/doc/write_coeff.txt b/doc/write_coeff.txt new file mode 100644 index 000000000..58624c9a8 --- /dev/null +++ b/doc/write_coeff.txt @@ -0,0 +1,46 @@ +"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c + +:link(lws,http://lammps.sandia.gov) +:link(ld,Manual.html) +:link(lc,Section_commands.html#comm) + +:line + +write_coeff command :h3 + +[Syntax:] + +write_coeff file :pre + +file = name of data file to write out :pre + +[Examples:] + +write_coeff polymer.coeff :pre + +[Description:] + +Write a text format file with the currently defined force field +coefficients in a way, that it can be read by LAMMPS with the +"include"_include.html command. In combination with the nocoeff +option of "write_data"_write_data.html this can be used to move +the Coeffs sections from a data file into a separate file. + +NOTE: The write_coeff command is not yet fully implemented in two +respects. First, some pair styles do not yet write their coefficient +information into the coeff file. This means you will need to specify +that information in your input script that reads the data file, via +the "pair_coeff"_pair_coeff.html command. + +:line + +[Restrictions:] + +none + +[Related commands:] + +"read_data"_read_data.html, "write_restart"_write_restart.html, +"write_data"_write_data.html + + diff --git a/src/write_coeff.cpp b/src/write_coeff.cpp new file mode 100644 index 000000000..9158e0c84 --- /dev/null +++ b/src/write_coeff.cpp @@ -0,0 +1,109 @@ +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + http://lammps.sandia.gov, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + 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. +------------------------------------------------------------------------- */ + +#include +#include +#include "write_coeff.h" +#include "pair.h" +#include "bond.h" +#include "angle.h" +#include "dihedral.h" +#include "improper.h" +#include "comm.h" +#include "force.h" +#include "universe.h" +#include "error.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- + called as write_coeff command in input script +------------------------------------------------------------------------- */ + +void WriteCoeff::command(int narg, char **arg) +{ + if (domain->box_exist == 0) + error->all(FLERR,"write_coeff command before simulation box is defined"); + + if (narg != 1) error->all(FLERR,"Illegal write_coeff command"); + + int n = strlen(arg[0]) + 5; + char *file = new char[n]; + + strcpy(file,"tmp."); + strcat(file,arg[0]); + + if (comm->me == 0) { + char str[256], coeff[256]; + FILE *one = fopen(file,"wb+"); + if (one == NULL) { + sprintf(str,"Cannot open coeff file %s",file); + error->one(FLERR,str); + } + + if (force->pair && force->pair->writedata) { + force->pair->init(); + fprintf(one,"# pair_style %s\npair_coeff\n",force->pair_style); + force->pair->write_data_all(one); + fprintf(one,"end\n"); + } + if (force->bond && force->bond->writedata) { + fprintf(one,"# bond_style %s\nbond_coeff\n",force->bond_style); + force->bond->write_data(one); + fprintf(one,"end\n"); + } + if (force->angle && force->angle->writedata) { + fprintf(one,"# angle_style %s\nangle_coeff\n",force->angle_style); + force->angle->write_data(one); + fprintf(one,"end\n"); + } + if (force->dihedral && force->dihedral->writedata) { + fprintf(one,"# dihedral_style %s\ndihedral_coeff\n",force->dihedral_style); + force->dihedral->write_data(one); + fprintf(one,"end\n"); + } + if (force->improper && force->improper->writedata) { + fprintf(one,"# improper_style %s\nimproper_coeff\n",force->improper_style); + force->improper->write_data(one); + fprintf(one,"end\n"); + } + rewind(one); + + FILE *two = fopen(file+4,"w"); + if (two == NULL) { + sprintf(str,"Cannot open coeff file %s",file+4); + error->one(FLERR,str); + } + fprintf(two,"# LAMMPS coeff file via write_coeff, version %s\n", + universe->version); + while(1) { + if (fgets(str,256,one) == NULL) break; + fputs(str,two); // style + fgets(str,256,one); // coeff + n = strlen(str); + strcpy(coeff,str); + coeff[n-1] = '\0'; + fgets(str,256,one); + while (strcmp(str,"end\n") != 0) { + fprintf(two,"%s %s",coeff,str); + fgets(str,256,one); + } + fputc('\n',two); + } + fclose(one); + fclose(two); + unlink(file); + } + + delete [] file; +} diff --git a/src/write_coeff.h b/src/write_coeff.h new file mode 100644 index 000000000..9b272b40f --- /dev/null +++ b/src/write_coeff.h @@ -0,0 +1,61 @@ +/* -*- c++ -*- ---------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + http://lammps.sandia.gov, Sandia National Laboratories + Steve Plimpton, sjplimp@sandia.gov + + 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. +------------------------------------------------------------------------- */ + +#ifdef COMMAND_CLASS + +CommandStyle(write_coeff,WriteCoeff) + +#else + +#ifndef LMP_WRITE_COEFF_H +#define LMP_WRITE_COEFF_H + +#include +#include "pointers.h" + +namespace LAMMPS_NS { + +class WriteCoeff : protected Pointers { + public: + WriteCoeff(class LAMMPS *lmp) : Pointers(lmp) {}; + void command(int, char **); +}; + +} + +#endif +#endif + +/* ERROR/WARNING messages: + +E: write_coeff command before simulation box is defined + +Self-explanatory. + +E: Illegal ... command + +Self-explanatory. Check the input script syntax and compare to the +documentation for the command. You can use -echo screen as a +command-line option when running LAMMPS to see the offending line. + +E: Atom count is inconsistent, cannot write data file + +The sum of atoms across processors does not equal the global number +of atoms. Probably some atoms have been lost. + +E: Cannot open coeff file %s + +The specified file cannot be opened. Check that the path and name are +correct. + +*/ diff --git a/test/02_commands/in.read_data-009 b/test/02_commands/in.read_data-009 new file mode 100644 index 000000000..c8257a61b --- /dev/null +++ b/test/02_commands/in.read_data-009 @@ -0,0 +1,16 @@ +# Solvated 5-mer peptide + +units real +atom_style full + +pair_style lj/charmm/coul/long 8.0 10.0 10.0 +bond_style harmonic +angle_style charmm +dihedral_style charmm +improper_style harmonic +kspace_style pppm 0.0001 + +read_data peptide.data + +write_coeff peptide.coeff +