diff --git a/src/atom.cpp b/src/atom.cpp index 2c2ebd911f..fb444e3e79 100644 --- a/src/atom.cpp +++ b/src/atom.cpp @@ -49,6 +49,7 @@ using namespace MathConst; #define DELTA 1 #define EPSILON 1.0e-6 +#define EPS_ZCOORD 1.0e-12 #define MAXLINE 256 /* ---------------------------------------------------------------------- @@ -1076,6 +1077,7 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset, // if periodic and I am lo/hi proc, adjust bounds by EPSILON // ensures all data atoms will be owned even with round-off + int dimension = domain->dimension; int triclinic = domain->triclinic; double epsilon[3]; @@ -1165,7 +1167,7 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset, imx = utils::inumeric(FLERR,values[iptr],false,lmp); imy = utils::inumeric(FLERR,values[iptr+1],false,lmp); imz = utils::inumeric(FLERR,values[iptr+2],false,lmp); - if ((domain->dimension == 2) && (imz != 0)) + if ((dimension == 2) && (imz != 0)) error->all(FLERR,"Z-direction image flag must be 0 for 2d-systems"); if ((!domain->xperiodic) && (imx != 0)) { reset_image_flag[0] = true; imx = 0; } if ((!domain->yperiodic) && (imy != 0)) { reset_image_flag[1] = true; imy = 0; } @@ -1179,6 +1181,15 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, tagint mol_offset, xdata[1] = utils::numeric(FLERR,values[xptr+1],false,lmp); xdata[2] = utils::numeric(FLERR,values[xptr+2],false,lmp); + // for 2d simulation, check if z coord is within EPS_ZCOORD of zero + // then set to zero + + if (dimension == 2) { + if (fabs(xdata[2]) > EPS_ZCOORD) + error->all(FLERR,"Read_data atom z coord is non-zero for 2d simulation"); + xdata[2] = 0.0; + } + // convert atom coords from general triclinic to restricted triclinic if (triclinic_general) domain->general_to_restricted_coords(xdata); diff --git a/src/create_atoms.cpp b/src/create_atoms.cpp index f7347f9ad1..0a8c462688 100644 --- a/src/create_atoms.cpp +++ b/src/create_atoms.cpp @@ -54,6 +54,8 @@ static constexpr double INV_P_CONST = 0.7548777; static constexpr double INV_SQ_P_CONST = 0.5698403; static constexpr int DEFAULT_MAXTRY = 1000; +#define EPS_ZCOORD 1.0e-12 + enum { BOX, REGION, SINGLE, RANDOM, MESH }; enum { ATOM, MOLECULE }; enum { COUNT, INSERT, INSERT_SELECTED }; @@ -1163,7 +1165,7 @@ void CreateAtoms::add_lattice() { // add atoms on general triclinic lattice if Domain has setting for it // verify lattice is valid for general triclinic - + int triclinic_general = domain->triclinic_general; if (triclinic_general) { @@ -1272,7 +1274,8 @@ void CreateAtoms::add_lattice() // decrement lo, increment hi to avoid round-off issues in lattice->bbox(), // which can lead to missing atoms in rare cases // extra decrement of lo if min < 0, since static_cast(-1.5) = -1 - + // for 2d simulation, klo = khi = 0 so just one plane of atoms + ilo = static_cast(xmin) - 1; jlo = static_cast(ymin) - 1; klo = static_cast(zmin) - 1; @@ -1284,7 +1287,7 @@ void CreateAtoms::add_lattice() if (ymin < 0.0) jlo--; if (zmin < 0.0) klo--; - printf("LOOP LATTICE bounds: %d %d: %d %d: %d %d\n",ilo,ihi,jlo,jhi,klo,khi); + if (domain->dimension == 2) klo = khi = 0; // count lattice sites on each proc @@ -1351,6 +1354,7 @@ void CreateAtoms::loop_lattice(int action) { int i, j, k, m; + int dimension = domain->dimension; int triclinic_general = domain->triclinic_general; const double *const *const basis = domain->lattice->basis; @@ -1372,9 +1376,18 @@ void CreateAtoms::loop_lattice(int action) domain->lattice->lattice2box(x[0], x[1], x[2]); // convert from general to restricted triclinic coords + // for 2d simulation, check if z coord is within EPS_ZCOORD of zero + // then set to zero + + if (triclinic_general) { + domain->general_to_restricted_coords(x); + if (dimension == 2) { + if (fabs(x[2]) > EPS_ZCOORD) + error->all(FLERR,"Create_atoms atom z coord is non-zero for 2d simulation"); + x[2] = 0.0; + } + } - if (triclinic_general) domain->general_to_restricted_coords(x); - // if a region was specified, test if atom is in it if (style == REGION) diff --git a/src/create_box.cpp b/src/create_box.cpp index 6a04744826..8cf6065962 100644 --- a/src/create_box.cpp +++ b/src/create_box.cpp @@ -63,6 +63,7 @@ void CreateBox::command(int narg, char **arg) int iarg = 2; if (region) { + // region is not prism // setup orthogonal box // set simulation domain from region extent @@ -94,6 +95,11 @@ void CreateBox::command(int narg, char **arg) domain->yz = prism->yz; } + if (domain->dimension == 2) { + if (domain->boxlo[2] >= 0.0 || domain->boxhi[2] <= 0.0) + error->all(FLERR,"Create_box region zlo/zhi for 2d simulation must straddle 0.0"); + } + // setup general triclinic box (with no region) // read next box extent arguments to create ABC edge vectors + origin // define_general_triclinic() converts @@ -114,7 +120,12 @@ void CreateBox::command(int narg, char **arg) double clo = utils::numeric(FLERR, arg[iarg + 4], false, lmp); double chi = utils::numeric(FLERR, arg[iarg + 5], false, lmp); iarg += 6; - + + if (domain->dimension == 2) { + if (clo >= 0.0 || clo <= 0.0) + error->all(FLERR,"Create_box region clo/chi for 2d simulation must straddle 0.0"); + } + // use lattice2box() to generate origin and ABC vectors // origin = abc lo // ABC vectors = hi in one dim - origin diff --git a/src/lattice.cpp b/src/lattice.cpp index 6d7450fd2b..ef2b4c344d 100644 --- a/src/lattice.cpp +++ b/src/lattice.cpp @@ -238,14 +238,14 @@ Lattice::Lattice(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp) if (dimension == 2) { if (origin[2] != 0.0) error->all(FLERR, - "Lattice settings are not compatible with 2d simulation"); + "Lattice origin not compatible with 2d simulation"); if (a1[2] != 0.0 || a2[2] != 0.0 || a3[0] != 0.0 || a3[1] != 0.0) error->all(FLERR, - "Lattice settings are not compatible with 2d simulation"); + "Lattice a1/a2/a3 vectors are not compatible with 2d simulation"); if (orientx[2] != 0 || orienty[2] != 0 || orientz[0] != 0 || orientz[1] != 0) error->all(FLERR, - "Lattice settings are not compatible with 2d simulation"); + "Lattice orient vectors are not compatible with 2d simulation"); } if (spaceflag) { diff --git a/src/read_data.cpp b/src/read_data.cpp index c19250c2aa..c59bec7476 100644 --- a/src/read_data.cpp +++ b/src/read_data.cpp @@ -491,11 +491,14 @@ void ReadData::command(int narg, char **arg) boxlo[0] = boxlo[1] = boxlo[2] = -0.5; boxhi[0] = boxhi[1] = boxhi[2] = 0.5; - avec[0] = avec[1] = avec[2] = 0.0; - bvec[0] = bvec[1] = bvec[2] = 0.0; - cvec[0] = cvec[1] = cvec[2] = 0.0; + xy = xz = yz = 0.0; + avec[0] = bvec[1] = cvec[2] = 1.0; + avec[1] = avec[2] = 0.0; + bvec[0] = bvec[2] = 0.0; + cvec[0] = cvec[1] = 0.0; abc_origin[0] = abc_origin[1] = abc_origin[2] = 0.0; + if (domain->dimension == 2) abc_origin[2] = -0.5; keyword[0] = '\0'; @@ -526,6 +529,19 @@ void ReadData::command(int narg, char **arg) error->all(FLERR,"Read_data header cannot specify simulation box lo/hi/tilt and ABC vectors"); triclinic = triclinic_general = 1; } + + // check if simulation box specified correctly for 2d + + if (domain->dimension == 2) { + if (triclinic_general == 0) { + if (boxlo[2] >= 0.0 || boxhi[2] <= 0.0) + error->all(FLERR,"Read_data zlo/zhi for 2d simulation must straddle 0.0"); + } else if (triclinic_general == 1) { + if (cvec[0] != 0.0 || cvec[1] != 0.0 || cvec[2] != 1.0 || abc_origin[2] != -0.5) + error->all(FLERR,"Read_data cvec and/or abc_origin is invalid for " + "2d simulation with general triclinic box"); + } + } // problem setup using info from header // only done once, if firstpass and first data file