ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

whitespace cleanup

Browse Source 
replaced all tabs with 8 spaces and deleted all trailing whitespace
This commit is contained in:
Axel Kohlmeyer
2016-03-03 16:19:23 -05:00
parent c144dd9765
commit 8bce407ab5
18 changed files with 112 additions and 112 deletions
Download Patch File Download Diff File Expand all files Collapse all files

6
src/USER-MANIFOLD/fix_manifoldforce.cpp
View File

@ -57,7 +57,7 @@ FixManifoldForce::FixManifoldForce(LAMMPS *lmp, int narg, char **arg) :
error->warning(FLERR,"Minimizing with fix manifoldforce without hftn or quickmin is fishy");
}
// Command is given as
// fix <name> <group> manifoldforce manifold_name manifold_args
if( narg < 5 ){
@ -65,7 +65,7 @@ FixManifoldForce::FixManifoldForce(LAMMPS *lmp, int narg, char **arg) :
}
const char *m_name = arg[3];
ptr_m = create_manifold(m_name,lmp,narg,arg);
// Construct manifold from factory:
if( !ptr_m ){
char msg[2048];
@ -165,7 +165,7 @@ void FixManifoldForce::post_force(int vflag)
f[i][0] -= dot*n[0] * invn2;
f[i][1] -= dot*n[1] * invn2;
f[i][2] -= dot*n[2] * invn2;
}
}
}

2
src/USER-MANIFOLD/fix_manifoldforce.h
View File

@ -56,7 +56,7 @@ class FixManifoldForce : public Fix {
void post_force(int);
void post_force_respa(int, int, int);
void min_post_force(int);
private:
user_manifold::manifold *ptr_m;

50
src/USER-MANIFOLD/fix_nve_manifold_rattle.cpp
View File

@ -83,7 +83,7 @@ FixNVEManifoldRattle::FixNVEManifoldRattle( LAMMPS *lmp, int &narg, char **arg,
{
if( lmp->citeme) lmp->citeme->add(cite_fix_nve_manifold_rattle);
if( narg < 6 ) error->all(FLERR, "Illegal fix nve/manifold/rattle command");
// Set all bits/settings:
time_integrate = 1;
dynamic_group_allow = 1;
@ -110,7 +110,7 @@ FixNVEManifoldRattle::FixNVEManifoldRattle( LAMMPS *lmp, int &narg, char **arg,
if( !tstrs || !tvars || !tstyle || !is_var ){
error->all(FLERR, "Error creating manifold arg arrays");
}
// Loop over manifold args:
for( int i = 0; i < nvars; ++i ){
int len = 0, offset = 0;
@ -136,7 +136,7 @@ FixNVEManifoldRattle::FixNVEManifoldRattle( LAMMPS *lmp, int &narg, char **arg,
}
ptr_m->post_param_init();
// Loop over rest of args:
int argi = 6 + nvars;
while( argi < narg ){
@ -144,14 +144,14 @@ FixNVEManifoldRattle::FixNVEManifoldRattle( LAMMPS *lmp, int &narg, char **arg,
nevery = force->inumeric(FLERR,arg[argi+1]);
argi += 2;
}else if( error_on_unknown_keyword ){
char msg[2048];
char msg[2048];
sprintf(msg,"Error parsing arg \"%s\".\n", arg[argi]);
error->all(FLERR, msg);
}else{
argi += 1;
}
}
}
@ -182,7 +182,7 @@ void FixNVEManifoldRattle::reset_dt()
dtf = 0.5 * update->dt * force->ftm2v;
}
void FixNVEManifoldRattle::print_stats( const char *header )
@ -241,7 +241,7 @@ int FixNVEManifoldRattle::setmask()
mask |= INITIAL_INTEGRATE;
mask |= FINAL_INTEGRATE;
if( nevery > 0 ) mask |= END_OF_STEP;
return mask;
}
@ -266,7 +266,7 @@ void FixNVEManifoldRattle::update_var_params()
stats.x_iters_per_atom = 0.0;
stats.v_iters_per_atom = 0.0;
}
double **ptr_params = ptr_m->get_params();
for( int i = 0; i < nvars; ++i ){
if( is_var[i] ){
@ -309,7 +309,7 @@ int FixNVEManifoldRattle::dof(int igroup)
// centre of mass corrections:
if( dofs <= 1 ) dofs = 0;
stats.dofs_removed = dofs;
return dofs;
}
@ -371,11 +371,11 @@ void FixNVEManifoldRattle::nve_x_rattle(int igroup, int groupbit)
int *mask = atom->mask;
int nlocal = atom->nlocal;
int natoms = 0;
if (igroup == atom->firstgroup){
nlocal = atom->nfirst;
}
if (rmass) {
for (int i = 0; i < nlocal; i++){
@ -470,7 +470,7 @@ void FixNVEManifoldRattle::rattle_manifold_x(double *x, double *v,
vo[0] = v[0];
vo[1] = v[1];
vo[2] = v[2];
xo[0] = x[0];
xo[1] = x[1];
xo[2] = x[2];
@ -490,10 +490,10 @@ void FixNVEManifoldRattle::rattle_manifold_x(double *x, double *v,
no_dt[2] = dtfm*no[2];
// Assume that no_nn is roughly constant during iteration:
const double c_inv = 1.0 / c;
while ( 1 ) {
v[0] = vt[0] - l*no_dt[0];
v[1] = vt[1] - l*no_dt[1];
@ -530,7 +530,7 @@ void FixNVEManifoldRattle::rattle_manifold_x(double *x, double *v,
++iters;
if( (res < tolerance) || (iters >= max_iter) ) break;
// Update nn and g.
gg = ptr_m->g(x);
ptr_m->n(x,nn);
@ -543,7 +543,7 @@ void FixNVEManifoldRattle::rattle_manifold_x(double *x, double *v,
tagi, x[0], x[1], x[2], res, iters);
error->one(FLERR,msg);
}
// "sync" x and v:
v[0] = vt[0] - l*no_dt[0];
v[1] = vt[1] - l*no_dt[1];
@ -564,7 +564,7 @@ void FixNVEManifoldRattle::rattle_manifold_v(double *v, double *f,
v[i][0] += dtfm * f[i][0];
v[i][1] += dtfm * f[i][1];
v[i][2] += dtfm * f[i][2];
Now you add the rattle-like update:
vold - vnew + dtfm * F + mu * n_new = 0
dot( vnew, n_new ) = 0
@ -604,7 +604,7 @@ void FixNVEManifoldRattle::rattle_manifold_v(double *v, double *f,
const double n_inv = 1.0 / nn2;
const double c_inv = 1.0 / c;
do{
R[0] = vt[0] - v[0] - l * no_dt[0];
R[1] = vt[1] - v[1] - l * no_dt[1];
@ -613,28 +613,28 @@ void FixNVEManifoldRattle::rattle_manifold_v(double *v, double *f,
// Analytic solution to system J*(dx,dy,dz,dl)^T = R
nn_R = n[0]*R[0] + n[1]*R[1] + n[2]*R[2];
dv[3] = -nn_R - R[3];
dv[3] *= n_inv;
dv[0] = -n[0]*dv[3] - R[0];
dv[1] = -n[1]*dv[3] - R[1];
dv[2] = -n[2]*dv[3] - R[2];
dv[3] *= c_inv;
v[0] -= dv[0];
v[1] -= dv[1];
v[2] -= dv[2];
l -= dv[3];
res = infnorm<4>(R);
++iters;
}while( (res > tolerance) && (iters < max_iter) );
if( iters >= max_iter && res >= tolerance ){
char msg[2048];
sprintf(msg,"Failed to constrain atom %d (x = (%f, %f, %f)! res = %e, iters = %d\n",
tagi, x[0], x[1], x[2], res, iters);
error->all(FLERR,msg);
char msg[2048];
sprintf(msg,"Failed to constrain atom %d (x = (%f, %f, %f)! res = %e, iters = %d\n",
tagi, x[0], x[1], x[2], res, iters);
error->all(FLERR,msg);
}
stats.v_iters += iters;

6
src/USER-MANIFOLD/fix_nve_manifold_rattle.h
View File

@ -64,7 +64,7 @@ namespace LAMMPS_NS {
int dofs_removed;
bigint last_out;
};
FixNVEManifoldRattle(LAMMPS *, int &, char **, int = 1);
virtual ~FixNVEManifoldRattle();
// All this stuff is interface, so you DO need to implement them.
@ -95,7 +95,7 @@ namespace LAMMPS_NS {
statistics stats;
int update_style;
int nvars;
user_manifold::manifold *ptr_m;
@ -124,7 +124,7 @@ 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: There is no manifold named ...
E: There is no manifold named ...
Self-explanatory. You requested a manifold whose name was not
registered at the factory.

36
src/USER-MANIFOLD/fix_nvt_manifold_rattle.cpp
View File

@ -90,10 +90,10 @@ FixNVTManifoldRattle::FixNVTManifoldRattle(LAMMPS *lmp, int narg, char **arg,
// Set all bits/settings:
dof_flag = 1;
dthalf = dt4 = dt8 = 0;
t_start = t_stop = t_period = t_current = t_target = ke_target = 0.0;
t_freq = drag = tdrag_factor = 0;
boltz = force->boltz, nktv2p = force->nktv2p;
tdof = 0;
mtchain = 3;
@ -117,11 +117,11 @@ FixNVTManifoldRattle::FixNVTManifoldRattle(LAMMPS *lmp, int narg, char **arg,
}else if( strcmp( arg[argi], "tchain" ) == 0 ){
if( argi+1 >= narg )
error->all(FLERR,"Keyword 'tchain' needs 1 argument");
mtchain = force->inumeric(FLERR, arg[argi+1]);
argi += 2;
}else if( error_on_unknown_keyword ){
char msg[2048];
char msg[2048];
sprintf(msg,"Error parsing arg \"%s\".\n", arg[argi]);
error->all(FLERR, msg);
}else{
@ -148,7 +148,7 @@ FixNVTManifoldRattle::FixNVTManifoldRattle(LAMMPS *lmp, int narg, char **arg,
newarg[1] = group->names[igroup];
newarg[2] = (char*) "temp";
modify->add_compute(3,newarg);
delete [] newarg;
int icompute = modify->find_compute(id_temp);
@ -159,7 +159,7 @@ FixNVTManifoldRattle::FixNVTManifoldRattle(LAMMPS *lmp, int narg, char **arg,
temperature = modify->compute[icompute];
if( temperature->tempbias ) which = BIAS;
else which = NOBIAS;
// Set t_freq from t_period
t_freq = 1.0 / t_period;
@ -174,8 +174,8 @@ FixNVTManifoldRattle::FixNVTManifoldRattle(LAMMPS *lmp, int narg, char **arg,
for( int ich = 0; ich < mtchain; ++ich ){
eta[ich] = eta_dot[ich] = eta_dotdot[ich] = 0.0;
}
}
/* ---------------------------------------------------------------------- */
@ -201,7 +201,7 @@ int FixNVTManifoldRattle::setmask()
mask |= INITIAL_INTEGRATE;
mask |= FINAL_INTEGRATE;
if( nevery > 0 ) mask |= END_OF_STEP;
return mask;
}
@ -241,9 +241,9 @@ void FixNVTManifoldRattle::setup(int vflag)
for( int ich = 1; ich < mtchain; ++ich ){
eta_mass[ich] = boltz * t_target * inv_t_freq2;
}
for( int ich = 1; ich < mtchain; ++ich ){
eta_dotdot[ich] = (eta_mass[ich-1]*eta_dot[ich-1]*eta_dot[ich-1] -
eta_dotdot[ich] = (eta_mass[ich-1]*eta_dot[ich-1]*eta_dot[ich-1] -
boltz * t_target ) / eta_mass[ich];
}
}
@ -252,7 +252,7 @@ void FixNVTManifoldRattle::compute_temp_target()
{
t_current = temperature->compute_scalar();
tdof = temperature->dof;
tdof = temperature->dof;
double delta = update->ntimestep - update->beginstep;
if (delta != 0.0){
@ -270,7 +270,7 @@ void FixNVTManifoldRattle::nhc_temp_integrate()
// t_current = temperature->compute_scalar();
// tdof = temperature->dof;
compute_temp_target();
double expfac, kecurrent = tdof * boltz * t_current;
double inv_t_freq2 = 1.0 / (t_freq*t_freq);
eta_mass[0] = tdof * boltz * t_target * inv_t_freq2;
@ -298,7 +298,7 @@ void FixNVTManifoldRattle::nhc_temp_integrate()
eta_dot[0] *= tdrag_factor * expfac;
factor_eta = exp(-dthalf*eta_dot[0]);
if( factor_eta == 0 ){
char msg[2048];
sprintf(msg, "WTF, factor_eta is 0! dthalf = %f, eta_dot[0] = %f",
@ -344,7 +344,7 @@ void FixNVTManifoldRattle::nh_v_temp()
if( which == NOBIAS ){
for( int i = 0; i < nlocal; ++i ){
if( mask[i] & groupbit ){
@ -363,7 +363,7 @@ void FixNVTManifoldRattle::nh_v_temp()
temperature->restore_bias(i,v[i]);
}
}
}
}
}
@ -394,7 +394,7 @@ void FixNVTManifoldRattle::final_integrate()
void FixNVTManifoldRattle::reset_dt()
{
FixNVEManifoldRattle::reset_dt();
dthalf = 0.5 * update->dt;
dt4 = 0.25 * update->dt;
dt8 = 0.125 * update->dt;
@ -410,6 +410,6 @@ double FixNVTManifoldRattle::memory_usage()
{
double bytes = FixNVEManifoldRattle::memory_usage();
bytes += (4*mtchain+1)*sizeof(double);
return bytes;
}

6
src/USER-MANIFOLD/fix_nvt_manifold_rattle.h
View File

@ -67,13 +67,13 @@ namespace LAMMPS_NS {
virtual double memory_usage();
protected:
void compute_temp_target();
void nhc_temp_integrate();
void nh_v_temp();
double dthalf, dt4, dt8;
char *id_temp;
@ -107,7 +107,7 @@ 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: There is no manifold named ...
E: There is no manifold named ...
Self-explanatory. You requested a manifold whose name was not
registered at the factory.

4
src/USER-MANIFOLD/manifold.h
View File

@ -71,7 +71,7 @@ namespace user_manifold {
protected:
double *params;
};
// Some utility functions that are templated, so I implement them
@ -90,7 +90,7 @@ namespace user_manifold {
inline double dot( double *a, double *b ){
return a[0]*b[0] + a[1]*b[1] + a[2]*b[2];
}
} // namespace LAMMPS_NS

4
src/USER-MANIFOLD/manifold_dumbbell.cpp
View File

@ -12,7 +12,7 @@ manifold_dumbbell::manifold_dumbbell( LAMMPS *lmp, int argc, char **argv ) : man
/*
* Equation for dumbbell is:
*
*
* -x^2 - y^2 + (a^2 - (z/c)^2)*( 1 + A*sin(const(z) *z^2) )^4
* const(z) = (z < 0) ? B2 : B
* params[4] = { a, A, B, c }
@ -37,7 +37,7 @@ void manifold_dumbbell::n( const double *x, double *nn )
nn[0] = -2.0*x[0];
nn[1] = -2.0*x[1];
nn[2] = 8.0*A*A*A*A*B*x[2]*(a*a-(x[2]*x[2]/(c*c)))*cos(B*x[2]*x[2]) *
pow(sin(B*x[2]*x[2]),3)-2*x[2]*(1+pow(A*sin(B*x[2]*x[2]),4))/(c*c);
pow(sin(B*x[2]*x[2]),3)-2*x[2]*(1+pow(A*sin(B*x[2]*x[2]),4))/(c*c);
}

4
src/USER-MANIFOLD/manifold_ellipsoid.h
View File

@ -18,8 +18,8 @@ namespace user_manifold {
static const char* type(){ return "ellipsoid"; }
virtual const char *id(){ return type(); }
static int expected_argc(){ return NPARAMS; }
virtual int nparams(){ return NPARAMS; }
static int expected_argc(){ return NPARAMS; }
virtual int nparams(){ return NPARAMS; }
};
}

4
src/USER-MANIFOLD/manifold_factory.cpp
View File

@ -65,8 +65,8 @@ manifold* LAMMPS_NS::user_manifold::create_manifold(const char *mname,
make_manifold_if<manifold_spine> ( &man, mname, lmp, narg, arg );
make_manifold_if<manifold_thylakoid> ( &man, mname, lmp, narg, arg );
make_manifold_if<manifold_torus> ( &man, mname, lmp, narg, arg );
return man;
}

2
src/USER-MANIFOLD/manifold_factory.h
View File

@ -101,7 +101,7 @@ namespace user_manifold {
manifold* create_manifold(const char *, LAMMPS *,
int , char ** );
} // namespace user_manifold
} // namespace LAMMPS_NS

8
src/USER-MANIFOLD/manifold_sphere.h
View File

@ -3,7 +3,7 @@
#include "manifold.h"
namespace LAMMPS_NS {
namespace LAMMPS_NS {
namespace user_manifold {
@ -29,7 +29,7 @@ namespace user_manifold {
nn[0] = 2*x[0];
nn[1] = 2*x[1];
nn[2] = 2*x[2];
return r2 - R*R;
}
@ -39,7 +39,7 @@ namespace user_manifold {
nn[1] = 2*x[1];
nn[2] = 2*x[2];
}
virtual void H( double *x, double h[3][3] )
{
h[0][1] = h[0][2] = h[1][0] = h[1][2] = h[2][0] = h[2][1] = 0.0;
@ -49,7 +49,7 @@ namespace user_manifold {
static const char* type(){ return "sphere"; }
virtual const char *id(){ return type(); }
static int expected_argc(){ return NPARAMS; }
virtual int nparams(){ return NPARAMS; }
virtual int nparams(){ return NPARAMS; }
};
}

4
src/USER-MANIFOLD/manifold_spine.cpp
View File

@ -11,7 +11,7 @@ manifold_spine::manifold_spine( LAMMPS *lmp, int argc, char **argv ) : manifold(
/*
* Equation for spine is:
*
*
* -x^2 - y^2 + (a^2 - (z/c)^2)*( 1 + A*sin(const(z) *z^2) )^4
* const(z) = (z < 0) ? B2 : B
* params[5] = { a, A, B, B2, c }
@ -43,7 +43,7 @@ void manifold_spine::n( const double *x, double *nn )
nn[0] = -2.0*x[0];
nn[1] = -2.0*x[1];
nn[2] = 8.0*A*A*A*A*B*x[2]*(a*a-(x[2]*x[2]/(c*c)))*cos(B*x[2]*x[2]) *
pow(sin(B*x[2]*x[2]),3)-2*x[2]*(1+pow(A*sin(B*x[2]*x[2]),4))/(c*c);
pow(sin(B*x[2]*x[2]),3)-2*x[2]*(1+pow(A*sin(B*x[2]*x[2]),4))/(c*c);
}else{
nn[0] = -2.0*x[0];
nn[1] = -2.0*x[1];

10
src/USER-MANIFOLD/manifold_supersphere.h
View File

@ -20,7 +20,7 @@ namespace user_manifold {
{
return (a > 0) - (a < 0);
}
virtual double g( const double *x )
{
double R = params[0];
@ -30,7 +30,7 @@ namespace user_manifold {
double zz = fabs(x[2]);
double rr = pow(xx,q) + pow(yy,q) + pow(zz,q);
return rr - pow(R,q);
}
@ -45,14 +45,14 @@ namespace user_manifold {
nn[1] = q * my_sign(x[1])*pow(yy,q-1);
nn[2] = q * my_sign(x[2])*pow(zz,q-1);
}
static const char* type(){ return "supersphere"; }
virtual const char *id(){ return type(); }
static int expected_argc(){ return NPARAMS; }
virtual int nparams(){ return NPARAMS; }
virtual int nparams(){ return NPARAMS; }
};
}
}
#endif
#endif

56
src/USER-MANIFOLD/manifold_thylakoid.cpp
View File

@ -119,7 +119,7 @@ void manifold_thylakoid::n( const double *x, double *n )
thyla_part *manifold_thylakoid::get_thyla_part( const double *x, int *err_flag, std::size_t *idx )
{
for( std::size_t i = 0; i < parts.size(); ++i ){
thyla_part *p = parts[i];
if( is_in_domain(p,x) ){
@ -146,8 +146,8 @@ void manifold_thylakoid::init_domains()
"wB = %f and lB = %f! %f > %f\n", LT, wB, lB, wB + 2*lB, LT);
error->one(FLERR,msg);
}
// Determine some constant coordinates:
// Determine some constant coordinates:
x0 = -( 0.5*LB + lB + lT + LT + lT + pad);
y0 = -( 0.5*LT + lT + pad );
z0 = -15;
@ -168,7 +168,7 @@ void manifold_thylakoid::init_domains()
char msg[2048];
sprintf(msg,"Thylakoid expects ylo of at most %f, but found %f",
y0, domain->boxlo[1]);
error->one(FLERR,msg);
error->one(FLERR,msg);
}
#endif
@ -178,7 +178,7 @@ void manifold_thylakoid::init_domains()
x1 = -x0;
y1 = -y0;
z1 = -z0;
Lx = x1 - x0;
Ly = y1 - y0;
Lz = z1 - z0;
@ -204,7 +204,7 @@ void manifold_thylakoid::init_domains()
// Create and add the manifold parts to the array.
thyla_part *p;
// Determine coordinates of domain boundaries and centres of "mass":
thyla_part_geom cllb, cllt, clrb, clrt; // Left thylakoid cylinder parts
@ -216,9 +216,9 @@ void manifold_thylakoid::init_domains()
thyla_part_geom bl, br, bc; // Bridge left, right connectors and cylinder.
// The bridge is three parts.
// 1. A connector between the right face of the left grana and a cylinder
// 2. A connector between the left face of the right grana and a cylinder
@ -250,7 +250,7 @@ void manifold_thylakoid::init_domains()
p = make_cyl_to_plane_part(bl.pt[0], rB, lB, -1, bl.pt);
set_domain(p, bl.lo, bl.hi);
parts.push_back(p);
// 2.:
br.pt[0] = 0.5*LB;
br.pt[1] = 0;
@ -273,9 +273,9 @@ void manifold_thylakoid::init_domains()
p = make_cyl_to_plane_part(br.pt[0], rB, lB, 1, br.pt);
set_domain(p, br.lo, br.hi);
parts.push_back(p);
// 3.:
// Cylinder in between:
bc.pt[0] = 0;
@ -285,7 +285,7 @@ void manifold_thylakoid::init_domains()
bc.lo[0] = bl.pt[0];
bc.lo[1] = -Reff;
bc.lo[2] = -Reff;
bc.hi[0] = br.pt[0];
bc.hi[1] = Reff;
bc.hi[2] = Reff;
@ -298,7 +298,7 @@ void manifold_thylakoid::init_domains()
bc.lo[0], bc.hi[0], bc.lo[1], bc.hi[1], bc.lo[2], bc.hi[2]);
}
#endif // MANIFOLD_THYLAKOID_DEBUG
parts.push_back(p);
@ -312,17 +312,17 @@ void manifold_thylakoid::init_domains()
cllb.hi[0] = cllb.pt[0];
cllb.hi[1] = cllb.pt[1];
cllb.hi[2] = z1;
p = make_cyl_part(1,1,0,cllb.pt,lT);
set_domain(p, cllb.lo, cllb.hi);
parts.push_back(p);
// left left top cylinder
cllt = cllb;
cllt.lo[1] = y1 - pad - lT;
cllt.hi[1] = y1;
cllt.pt[1] = cllb.pt[1] + LT;
p = make_cyl_part(1,1,0,cllt.pt,lT);
set_domain(p, cllt.lo, cllt.hi);
parts.push_back(p);
@ -342,7 +342,7 @@ void manifold_thylakoid::init_domains()
clrt.pt[1] += LT;
clrt.lo[1] = y1 - pad - lT;
clrt.hi[1] = y1;
p = make_cyl_part(1,1,0,clrt.pt,lT);
set_domain(p, clrt.lo, clrt.hi);
parts.push_back(p);
@ -362,7 +362,7 @@ void manifold_thylakoid::init_domains()
p = make_plane_part(1,0,0,pll.pt);
set_domain(p, pll.lo, pll.hi);
parts.push_back(p);
// left bottom plane
plb.pt[0] = x0 + pad + lT + 0.5*LT;
plb.pt[1] = y0 + pad;
@ -373,17 +373,17 @@ void manifold_thylakoid::init_domains()
plb.hi[0] = plb.lo[0] + LT;
plb.hi[1] = plb.lo[1] + pad + lT;
plb.hi[2] = z1;
p = make_plane_part(0,1,0,plb.pt);
set_domain(p, plb.lo, plb.hi);
parts.push_back(p);
// left top plane
plt = plb;
plt.lo[1] = cllb.pt[1] + LT;
plt.hi[1] = y1;
plt.pt[1] = y1 - pad;
p = make_plane_part(0,1,0,plt.pt);
set_domain(p, plt.lo, plt.hi);
parts.push_back(p);
@ -436,7 +436,7 @@ void manifold_thylakoid::init_domains()
parts.push_back(p);
thyla_part_geom::mirror(thyla_part_geom::DIR_X, &prl , &pll );
p = make_plane_part(1,0,0,prl.pt);
set_domain(p, prl.lo, prl.hi);
@ -503,7 +503,7 @@ void manifold_thylakoid::set_domain( thyla_part *p, const std::vector<double> &l
p->xlo = lo[0];
p->ylo = lo[1];
p->zlo = lo[2];
p->xhi = hi[0];
p->yhi = hi[1];
p->zhi = hi[2];
@ -512,13 +512,13 @@ void manifold_thylakoid::set_domain( thyla_part *p, const std::vector<double> &l
int manifold_thylakoid::is_in_domain( thyla_part *part, const double *x )
{
bool domain_ok = (x[0] >= part->xlo) && (x[0] <= part->xhi) &&
(x[1] >= part->ylo) && (x[1] <= part->yhi) &&
(x[2] >= part->zlo) && (x[2] <= part->zhi);
(x[1] >= part->ylo) && (x[1] <= part->yhi) &&
(x[2] >= part->zlo) && (x[2] <= part->zhi);
if( !domain_ok ) return false;
// From here on out, domain is ok.
if( part->type == thyla_part::THYLA_TYPE_CYL_TO_PLANE ){
double R0 = part->params[1];
@ -528,8 +528,8 @@ int manifold_thylakoid::is_in_domain( thyla_part *part, const double *x )
double dist2 = y*y + z*z;
double RR = R0+R;
double RR2 = RR*RR;
if( dist2 < RR2 ){
return true;
}else{

6
src/USER-MANIFOLD/manifold_thylakoid.h
View File

@ -17,16 +17,16 @@ namespace user_manifold {
enum { NPARAMS = 3 };
manifold_thylakoid( LAMMPS *lmp, int, char ** );
virtual ~manifold_thylakoid();
virtual double g( const double *x );
virtual void n( const double *x, double *n );
static const char* type(){ return "thylakoid"; }
virtual const char *id(){ return type(); }
static int expected_argc(){ return NPARAMS; }
static int expected_argc(){ return NPARAMS; }
virtual int nparams(){ return NPARAMS; }
virtual void post_param_init();
virtual void checkup(); // Some diagnostics...
private:

12
src/USER-MANIFOLD/manifold_thylakoid_shared.cpp
View File

@ -60,7 +60,7 @@ thyla_part::thyla_part( int type, double *args, double xlo, double ylo, double z
*
* As params it takes X0, R0, R, and a sign that determines whether it is
* "left" (args[3] < 0) or "right" (args[3] > 0).
*
*
* The args are: X0, R0, R, x0, y0, z0, sign.
*/
params[0] = args[0];
@ -114,7 +114,7 @@ double thyla_part::g(const double *x)
double dx = x[0] - X0;
double dy = x[1] - Y0;
double dz = x[2] - Z0;
return a*dx*dx + b*dy*dy + c*dz*dz - R*R;
break;
}
@ -122,7 +122,7 @@ double thyla_part::g(const double *x)
double X0 = params[0];
double R0 = params[1];
double R = params[2];
// Determine the centre of the sphere.
double dx = (x[0] - X0);
double dyz = sqrt( x[1]*x[1] + x[2]*x[2] );
@ -184,16 +184,16 @@ void thyla_part::n( const double *x, double *n )
double R0 = params[1];
double R = params[2];
double s = (params[6] > 0.0) ? 1.0 : -1.0;
// Determine the centre of the sphere.
double dx = s*(x[0] - X0);
double ryz = sqrt( x[1]*x[1] + x[2]*x[2] );
// Maybe sign is important here...
// Normalize g and n so that the normal is continuous:
double norm = 1.0 / (2.0 * R);
n[0] = s*2*dx*norm;
double const_part = 1.0 - (R0 + R) / ryz;
n[1] = 2*x[1]*const_part*norm;
n[2] = 2*x[2]*const_part*norm;

4
src/USER-MANIFOLD/manifold_thylakoid_shared.h
View File

@ -34,7 +34,7 @@ namespace user_manifold {
int err_flag;
tagint a_id;
double xlo, xhi, ylo, yhi, zlo, zhi;
double x0, y0, z0;
@ -50,7 +50,7 @@ namespace user_manifold {
enum DIRS { DIR_X, DIR_Y, DIR_Z };
static void mirror( unsigned int axis, thyla_part_geom *m,
const thyla_part_geom *o );
}; // struct thyla_part_geom