diff --git a/lib/atc/ATC_Coupling.cpp b/lib/atc/ATC_Coupling.cpp
index c6e239cc40130931f5d92ee4ffb2f415edde6dfe..c5df4d741c3eec909719c2e3d3adb89be136ada7 100644
--- a/lib/atc/ATC_Coupling.cpp
+++ b/lib/atc/ATC_Coupling.cpp
@@ -1070,6 +1070,16 @@ namespace ATC {
rhs);
}
}
+
+ //--------------------------------------------------
+ bool ATC_Coupling::reset_methods() const
+ {
+ bool resetMethods = ATC_Method::reset_methods() || atomicRegulator_->need_reset();
+ for (_ctiIt_ = timeIntegrators_.begin(); _ctiIt_ != timeIntegrators_.end(); ++_ctiIt_) {
+ resetMethods |= (_ctiIt_->second)->need_reset();
+ }
+ return resetMethods;
+ }
//--------------------------------------------------
void ATC_Coupling::initialize()
{
@@ -1122,10 +1132,120 @@ namespace ATC {
}
}
}
-
// prepare computes for first timestep
lammpsInterface_->computes_addstep(lammpsInterface_->ntimestep()+1);
+
+ // resetting precedence:
+ // time integrator -> kinetostat/thermostat -> time filter
+ // init_filter uses fieldRateNdFiltered which comes from the time integrator,
+ // which is why the time integrator is initialized first
+
+ // other initializations
+ if (reset_methods()) {
+ for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
+ (_tiIt_->second)->initialize();
+ }
+ atomicRegulator_->initialize();
+ }
+ extrinsicModelManager_.initialize();
+ if (timeFilterManager_.need_reset()) {// reset thermostat power
+ init_filter();
+ }
+ // clears need for reset
+ timeFilterManager_.initialize();
+ ghostManager_.initialize();
+
+ if (!initialized_) {
+ // initialize sources based on initial FE temperature
+ double dt = lammpsInterface_->dt();
+ prescribedDataMgr_->set_sources(time()+0.5*dt,sources_);
+ extrinsicModelManager_.set_sources(fields_,extrinsicSources_);
+ atomicRegulator_->compute_boundary_flux(fields_);
+ compute_atomic_sources(fieldMask_,fields_,atomicSources_);
+
+ // read in field data if necessary
+ if (useRestart_) {
+ RESTART_LIST data;
+ read_restart_data(restartFileName_,data);
+ useRestart_ = false;
+ }
+
+ // set consistent initial conditions, if requested
+ if (!timeFilterManager_.filter_dynamics() && consistentInitialization_) {
+
+ const INT_ARRAY & nodeType(nodalGeometryType_->quantity());
+
+ if (fieldSizes_.find(VELOCITY) != fieldSizes_.end()) {
+ DENS_MAT & velocity(fields_[VELOCITY].set_quantity());
+ DENS_MAN * nodalAtomicVelocity(interscaleManager_.dense_matrix("NodalAtomicVelocity"));
+ const DENS_MAT & atomicVelocity(nodalAtomicVelocity->quantity());
+ for (int i = 0; i<nNodes_; ++i) {
+
+ if (nodeType(i,0)==MD_ONLY) {
+ for (int j = 0; j < nsd_; j++) {
+ velocity(i,j) = atomicVelocity(i,j);
+ }
+ }
+ }
+ }
+
+ if (fieldSizes_.find(TEMPERATURE) != fieldSizes_.end()) {
+ DENS_MAT & temperature(fields_[TEMPERATURE].set_quantity());
+ DENS_MAN * nodalAtomicTemperature(interscaleManager_.dense_matrix("NodalAtomicTemperature"));
+ const DENS_MAT & atomicTemperature(nodalAtomicTemperature->quantity());
+
+ for (int i = 0; i<nNodes_; ++i) {
+
+ if (nodeType(i,0)==MD_ONLY) {
+ temperature(i,0) = atomicTemperature(i,0);
+ }
+ }
+ }
+
+ if (fieldSizes_.find(DISPLACEMENT) != fieldSizes_.end()) {
+ DENS_MAT & displacement(fields_[DISPLACEMENT].set_quantity());
+ DENS_MAN * nodalAtomicDisplacement(interscaleManager_.dense_matrix("NodalAtomicDisplacement"));
+ const DENS_MAT & atomicDisplacement(nodalAtomicDisplacement->quantity());
+ for (int i = 0; i<nNodes_; ++i) {
+
+ if (nodeType(i,0)==MD_ONLY) {
+ for (int j = 0; j < nsd_; j++) {
+ displacement(i,j) = atomicDisplacement(i,j);
+ }
+ }
+ }
+ }
+
+ //WIP_JAT update next two when full species time integrator is added
+ if (fieldSizes_.find(MASS_DENSITY) != fieldSizes_.end()) {
+ DENS_MAT & massDensity(fields_[MASS_DENSITY].set_quantity());
+ const DENS_MAT & atomicMassDensity(atomicFields_[MASS_DENSITY]->quantity());
+ for (int i = 0; i<nNodes_; ++i) {
+
+ if (nodeType(i,0)==MD_ONLY) {
+ massDensity(i,0) = atomicMassDensity(i,0);
+ }
+ }
+ }
+
+ if (fieldSizes_.find(SPECIES_CONCENTRATION) != fieldSizes_.end()) {
+ DENS_MAT & speciesConcentration(fields_[SPECIES_CONCENTRATION].set_quantity());
+ const DENS_MAT & atomicSpeciesConcentration(atomicFields_[SPECIES_CONCENTRATION]->quantity());
+ for (int i = 0; i<nNodes_; ++i) {
+
+ if (nodeType(i,0)==MD_ONLY) {
+ for (int j = 0; j < atomicSpeciesConcentration.nCols(); ++j) {
+ speciesConcentration(i,j) = atomicSpeciesConcentration(i,j);
+ }
+ }
+ }
+ }
+ }
+
+ initialized_ = true;
+ }
+
}
//-------------------------------------------------------------------
void ATC_Coupling::construct_time_integration_data()
@@ -1731,22 +1851,23 @@ namespace ATC {
}
//--------------------------------------------------------
- // mid_init_integrate
- // time integration between the velocity update and
- // the position lammps update of Verlet step 1
+ // init_integrate
+ // time integration of lammps atomic quantities
//--------------------------------------------------------
- void ATC_Coupling::mid_init_integrate()
+ void ATC_Coupling::init_integrate()
{
- double dt = lammpsInterface_->dt();
-
- // Compute nodal velocity at n+1/2, if needed
- for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
- (_tiIt_->second)->mid_initial_integrate1(dt);
- }
-
- atomicRegulator_->apply_mid_predictor(dt,lammpsInterface_->ntimestep());
-
- extrinsicModelManager_.mid_init_integrate();
+ atomTimeIntegrator_->init_integrate_velocity(dt());
+ ghostManager_.init_integrate_velocity(dt());
+ // account for other fixes doing time integration
+ interscaleManager_.fundamental_force_reset(LammpsInterface::ATOM_VELOCITY);
+
+ // apply constraints to velocity
+ atomicRegulator_->apply_mid_predictor(dt(),lammpsInterface_->ntimestep());
+
+ atomTimeIntegrator_->init_integrate_position(dt());
+ ghostManager_.init_integrate_position(dt());
+ // account for other fixes doing time integration
+ interscaleManager_.fundamental_force_reset(LammpsInterface::ATOM_POSITION);
}
///--------------------------------------------------------
@@ -1941,6 +2062,7 @@ namespace ATC {
output();
lammpsInterface_->computes_addstep(lammpsInterface_->ntimestep()+1); // adds next step to computes
+ //ATC_Method::post_final_integrate();
}
//=================================================================
diff --git a/lib/atc/ATC_Coupling.h b/lib/atc/ATC_Coupling.h
index 7b9ff149836124d29628369933421471e0ee779b..d9ac1485428cc3703ca0a96332ecd291ee502bf6 100644
--- a/lib/atc/ATC_Coupling.h
+++ b/lib/atc/ATC_Coupling.h
@@ -61,27 +61,25 @@ namespace ATC {
/** pre integration run */
virtual void initialize();
+ /** flags whether a methods reset is required */
+ virtual bool reset_methods() const;
+
/** post integration run : called at end of run or simulation */
virtual void finish();
/** first time, before atomic integration */
virtual void pre_init_integrate();
- /** Predictor phase, executed between velocity and position Verlet */
- virtual void mid_init_integrate();
+ /** Predictor phase, Verlet first step for velocity and position */
+ virtual void init_integrate();
/** Predictor phase, executed after Verlet */
virtual void post_init_integrate();
-#ifdef OBSOLETE
- /** Corrector phase, executed before Verlet */
- virtual void pre_final_integrate(){};
-#endif
+
/** Corrector phase, executed after Verlet*/
virtual void post_final_integrate();
-#ifdef OBSOLETE
-{lammpsInterface_->computes_addstep(lammpsInterface_->ntimestep()+1);};
-#endif
+
/** pre/post atomic force calculation in minimize */
virtual void min_pre_force(){};
virtual void min_post_force(){};
@@ -446,6 +444,8 @@ namespace ATC {
FIELDS rhs_; // for pde
FIELDS rhsAtomDomain_; // for thermostat
FIELDS boundaryFlux_; // for thermostat & rhs pde
+ // DATA structures for tracking individual species and molecules
+ FIELD_POINTERS atomicFields_;
/*@}*/
// workspace variables
diff --git a/lib/atc/ATC_CouplingEnergy.cpp b/lib/atc/ATC_CouplingEnergy.cpp
index 4592614782795de9ec7a2acec403d9d695d8ac1c..8ba10c4f8a08cbbbdbcfab461e3a1f2282ca479d 100644
--- a/lib/atc/ATC_CouplingEnergy.cpp
+++ b/lib/atc/ATC_CouplingEnergy.cpp
@@ -30,9 +30,6 @@ namespace ATC {
string matParamFile,
ExtrinsicModelType extrinsicModel)
: ATC_Coupling(groupName,perAtomArray,thisFix),
-#ifdef OBSOLETE
- nodalAtomicHeatCapacity_(NULL),
-#endif
nodalAtomicKineticTemperature_(NULL),
nodalAtomicConfigurationalTemperature_(NULL)
{
@@ -69,9 +66,6 @@ namespace ATC {
extVector_ = 1;
if (extrinsicModel != NO_MODEL)
sizeVector_ += extrinsicModelManager_.size_vector(sizeVector_);
-
- // create PE per atom ccompute
- //lammpsInterface_->create_compute_pe_peratom();
}
//--------------------------------------------------------
@@ -91,60 +85,6 @@ namespace ATC {
{
// Base class initalizations
ATC_Coupling::initialize();
-
- // resetting precedence:
- // time integrator -> thermostat -> time filter
- // init_filter uses fieldRateNdFiltered which comes from the time integrator,
- // which is why the time integrator is initialized first
-
- // other initializations
- if (reset_methods()) {
- for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
- (_tiIt_->second)->initialize();
- }
- atomicRegulator_->initialize();
- }
- extrinsicModelManager_.initialize();
- // reset thermostat power for time filter initial conditions for special cases
- if (timeFilterManager_.need_reset()) {
- init_filter();
- }
- // clears need for reset
- timeFilterManager_.initialize();
- ghostManager_.initialize();
-
- if (!initialized_) {
- // initialize sources based on initial FE temperature
- double dt = lammpsInterface_->dt();
- prescribedDataMgr_->set_sources(time()+0.5*dt,sources_);
- extrinsicModelManager_.set_sources(fields_,extrinsicSources_);
- atomicRegulator_->compute_boundary_flux(fields_);
- compute_atomic_sources(fieldMask_,fields_,atomicSources_);
-
- // read in field data if necessary
- if (useRestart_) {
- RESTART_LIST data;
- read_restart_data(restartFileName_,data);
- useRestart_ = false;
- }
-
- // set consistent initial conditions, if requested
- if (!timeFilterManager_.filter_dynamics()) {
- if (consistentInitialization_) {
-
- DENS_MAT & temperature(fields_[TEMPERATURE].set_quantity());
- DENS_MAN * nodalAtomicTemperature(interscaleManager_.dense_matrix("NodalAtomicTemperature"));
- const DENS_MAT & atomicTemperature(nodalAtomicTemperature->quantity());
- const INT_ARRAY & nodeType(nodalGeometryType_->quantity());
- for (int i = 0; i<nNodes_; ++i) {
-
- if (nodeType(i,0)==MD_ONLY)
- temperature(i,0) = atomicTemperature(i,0);
- }
- }
- }
- initialized_ = true;
- }
// reset integration field mask
intrinsicMask_.reset(NUM_FIELDS,NUM_FLUX);
@@ -257,21 +197,7 @@ namespace ATC {
TEMPERATURE);
interscaleManager_.add_dense_matrix(nodalAtomicTemperature,
"NodalAtomicTemperature");
-#ifdef OBSOLETE
- if (!useFeMdMassMatrix_) {
- // classical thermodynamic heat capacity of the atoms
- HeatCapacity * heatCapacity = new HeatCapacity(this);
- interscaleManager_.add_per_atom_quantity(heatCapacity,
- "AtomicHeatCapacity");
-
- // atomic thermal mass matrix
- nodalAtomicHeatCapacity_ = new AtfShapeFunctionRestriction(this,
- heatCapacity,
- shpFcn_);
- interscaleManager_.add_dense_matrix(nodalAtomicHeatCapacity_,
- "NodalAtomicHeatCapacity");
- }
-#endif
+
for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
(_tiIt_->second)->construct_transfers();
}
@@ -314,20 +240,7 @@ namespace ATC {
}
}
}
-#ifdef OBSOLETE
- //---------------------------------------------------------
- // compute_md_mass_matrix
- // compute the mass matrix arising from only atomistic
- // quadrature and contributions as a summation
- //---------------------------------------------------------
- void ATC_CouplingEnergy::compute_md_mass_matrix(FieldName thisField,
- DIAG_MAT & massMat)
- {
-
- if (thisField == TEMPERATURE)
- massMat.reset(nodalAtomicHeatCapacity_->quantity());
- }
-#endif
+
//--------------------------------------------------------
// modify
// parses inputs and modifies state of the filter
diff --git a/lib/atc/ATC_CouplingEnergy.h b/lib/atc/ATC_CouplingEnergy.h
index 7ce6063dced3765c5cd364de2be50ed46067f5e5..fce1a4550171bd3132a39e910533d50f121b5263 100644
--- a/lib/atc/ATC_CouplingEnergy.h
+++ b/lib/atc/ATC_CouplingEnergy.h
@@ -42,15 +42,6 @@ namespace ATC {
/** pre time integration */
virtual void initialize();
- /** flags whether a methods reset is required */
-
- virtual bool reset_methods() const {
- bool resetMethods = ATC_Method::reset_methods() || atomicRegulator_->need_reset();
- _ctiIt_ = timeIntegrators_.find(TEMPERATURE);
- if (_ctiIt_ == timeIntegrators_.end()) return resetMethods;
- return resetMethods || (_ctiIt_->second)->need_reset();
- };
-
/** compute vector for output */
virtual double compute_vector(int n);
@@ -69,14 +60,7 @@ namespace ATC {
/** sets the position/velocity of the ghost atoms */
virtual void set_ghost_atoms(){};
-#ifdef OBSOLETE
- /** compute the mass matrix components coming from MD integration */
- virtual void compute_md_mass_matrix(FieldName thisField,
- DIAG_MAT & massMats);
- /** operator to compute mass matrix from MD */
- AtfShapeFunctionRestriction * nodalAtomicHeatCapacity_;
-#endif
/** physics specific filter initialization */
void init_filter();
diff --git a/lib/atc/ATC_CouplingMass.cpp b/lib/atc/ATC_CouplingMass.cpp
index 7bc7e7e3c8281428f6c2c87e4b9574d40689daa7..801752d569082aa9ec4973e1183911dc3e490f24 100644
--- a/lib/atc/ATC_CouplingMass.cpp
+++ b/lib/atc/ATC_CouplingMass.cpp
@@ -123,68 +123,6 @@ namespace ATC {
// Base class initalizations
ATC_Coupling::initialize();
- // check that only all atoms
-
- if (bndyIntType_ != NO_QUADRATURE) throw ATC_Error("ATC_CouplingMass: only all atoms simulations are supported");
-
- // set consistent initial conditions, if requested
- if (!timeFilterManager_.filter_dynamics()) {
- if (consistentInitialization_) {
-
- DENS_MAT & massDensity(fields_[MASS_DENSITY].set_quantity());
- const DENS_MAT & atomicMassDensity(atomicFields_[MASS_DENSITY]->quantity());
-
- DENS_MAT & speciesConcentration(fields_[SPECIES_CONCENTRATION].set_quantity());
- const DENS_MAT & atomicSpeciesConcentration(atomicFields_[SPECIES_CONCENTRATION]->quantity());
-
- const INT_ARRAY & nodeType(nodalGeometryType_->quantity());
- for (int i = 0; i<nNodes_; ++i) {
-
- if (nodeType(i,0)==MD_ONLY) {
- massDensity(i,0) = atomicMassDensity(i,0);
- for (int j = 0; j < atomicSpeciesConcentration.nCols(); ++j) {
- speciesConcentration(i,j) = atomicSpeciesConcentration(i,j);
- }
- }
- }
- }
- }
-
-
- // other initializatifields_[SPECIES_CONCENTRATION].quantity()ons
- if (reset_methods()) {
- for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
- (_tiIt_->second)->initialize();
- }
- }
- extrinsicModelManager_.initialize(); // always needed to construct new Poisson solver
- if (timeFilterManager_.need_reset()) {
- init_filter();
- }
- // clears need for reset
- timeFilterManager_.initialize();
- atomicRegulator_->initialize();
- ghostManager_.initialize();
-
- if (!initialized_) {
- // initialize sources based on initial FE temperature
- double dt = lammpsInterface_->dt();
- // set sources
-
- prescribedDataMgr_->set_sources(time()+0.5*dt,sources_);
- extrinsicModelManager_.set_sources(fields_,extrinsicSources_);
- compute_atomic_sources(fieldMask_,fields_,atomicSources_);
-
- // read in field data if necessary
- if (useRestart_) {
- RESTART_LIST data;
- read_restart_data(restartFileName_,data);
- useRestart_ = false;
- }
-
- initialized_ = true;
- }
-
// reset integration field mask
intrinsicMask_.reset(NUM_FIELDS,NUM_FLUX);
intrinsicMask_ = false;
@@ -245,17 +183,7 @@ namespace ATC {
ATC_Coupling::init_filter();
}
-#ifdef OBSOLETE
- void ATC_CouplingMass::compute_md_mass_matrix(FieldName thisField,
- DIAG_MAT & massMat)
- {
- if (thisField == MASS_DENSITY ||
- thisField == SPECIES_CONCENTRATION) {
- massMat.reset(nodalAtomicVolume_->quantity());
- }
- }
-#endif
//WIP_JAT consolidate to coupling when we handle the temperature correctly
//--------------------------------------------------------
// pre_exchange
diff --git a/lib/atc/ATC_CouplingMass.h b/lib/atc/ATC_CouplingMass.h
index 582845a8a2eb415c8c28dd4d1132854a9e312527..1589a7a655a1780076f2209bda243f86f4fbed47 100644
--- a/lib/atc/ATC_CouplingMass.h
+++ b/lib/atc/ATC_CouplingMass.h
@@ -72,17 +72,10 @@ namespace ATC {
/** sets the position/velocity of the ghost atoms */
virtual void set_ghost_atoms(){};
-#ifdef OBSOLETE
- /** compute the mass matrix components coming from MD integration */
- virtual void compute_md_mass_matrix(FieldName thisField,
- DIAG_MAT & massMats);
-#endif
+
/** physics specific filter initialization */
void init_filter();
- // DATA structures for tracking individual species and molecules
- FIELD_POINTERS atomicFields_;
-
bool resetNlocal_;
diff --git a/lib/atc/ATC_CouplingMomentum.cpp b/lib/atc/ATC_CouplingMomentum.cpp
index 6eb8529a75f3cb5ebc8fedb703e9558b20516495..f59aa525cbfb06512d72aed78a65f0b032415278 100644
--- a/lib/atc/ATC_CouplingMomentum.cpp
+++ b/lib/atc/ATC_CouplingMomentum.cpp
@@ -33,10 +33,6 @@ namespace ATC {
PhysicsType intrinsicModel,
ExtrinsicModelType extrinsicModel)
: ATC_Coupling(groupName,perAtomArray,thisFix),
-#ifdef OBSOLETE
- nodalAtomicMass_(NULL),
- nodalAtomicCount_(NULL),
-#endif
refPE_(0)
{
// Allocate PhysicsModel
@@ -110,73 +106,6 @@ namespace ATC {
// Base class initalizations
ATC_Coupling::initialize();
- // check resetting precedence:
- // time integrator -> kinetostat -> time filter
-
- // other initializations
- if (reset_methods()) {
- for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
- (_tiIt_->second)->initialize();
- }
- atomicRegulator_->initialize();
- }
- extrinsicModelManager_.initialize();
- if (timeFilterManager_.need_reset()) { // reset kinetostat power
- init_filter();
- }
- // clears need for reset
- timeFilterManager_.initialize();
- ghostManager_.initialize();
-
- if (!initialized_) {
- // initialize sources based on initial FE temperature
- double dt = lammpsInterface_->dt();
- prescribedDataMgr_->set_sources(time()+0.5*dt,sources_);
- extrinsicModelManager_.set_sources(fields_,extrinsicSources_);
- atomicRegulator_->compute_boundary_flux(fields_);
- compute_atomic_sources(fieldMask_,fields_,atomicSources_);
-
- // read in field data if necessary
- if (useRestart_) {
- RESTART_LIST data;
- read_restart_data(restartFileName_,data);
- useRestart_ = false;
- }
-
- // set consistent initial conditions, if requested
- if (!timeFilterManager_.filter_dynamics()) {
- if (consistentInitialization_) {
-
- DENS_MAT & velocity(fields_[VELOCITY].set_quantity());
- DENS_MAN * nodalAtomicVelocity(interscaleManager_.dense_matrix("NodalAtomicVelocity"));
- const DENS_MAT & atomicVelocity(nodalAtomicVelocity->quantity());
- const INT_ARRAY & nodeType(nodalGeometryType_->quantity());
- for (int i = 0; i<nNodes_; ++i) {
-
- if (nodeType(i,0)==MD_ONLY) {
- for (int j = 0; j < nsd_; j++) {
- velocity(i,j) = atomicVelocity(i,j);
- }
- }
- }
- if (trackDisplacement_) {
- DENS_MAT & displacement(fields_[DISPLACEMENT].set_quantity());
- DENS_MAN * nodalAtomicDisplacement(interscaleManager_.dense_matrix("NodalAtomicDisplacement"));
- const DENS_MAT & atomicDisplacement(nodalAtomicDisplacement->quantity());
- for (int i = 0; i<nNodes_; ++i) {
-
- if (nodeType(i,0)==MD_ONLY) {
- for (int j = 0; j < nsd_; j++) {
- displacement(i,j) = atomicDisplacement(i,j);
- }
- }
- }
- }
- }
- }
- initialized_ = true;
- }
-
// reset integration field mask
intrinsicMask_.reset(NUM_FIELDS,NUM_FLUX);
intrinsicMask_ = false;
@@ -252,27 +181,7 @@ namespace ATC {
interscaleManager_.add_dense_matrix(nodalAtomicDisplacement,
"NodalAtomicDisplacement");
}
-#ifdef OBSOLETE
- // atomic mass matrix data
- if (!useFeMdMassMatrix_) {
- // atomic momentum mass matrix
- FundamentalAtomQuantity * atomicMass = interscaleManager_.fundamental_atom_quantity(LammpsInterface::ATOM_MASS);
- nodalAtomicMass_ = new AtfShapeFunctionRestriction(this,
- atomicMass,
- shpFcn_);
- interscaleManager_.add_dense_matrix(nodalAtomicMass_,
- "AtomicMomentumMassMat");
-
- // atomic dimensionless mass matrix
- ConstantQuantity<double> * atomicOnes = new ConstantQuantity<double>(this,1);
- interscaleManager_.add_per_atom_quantity(atomicOnes,"AtomicOnes");
- nodalAtomicCount_ = new AtfShapeFunctionRestriction(this,
- atomicOnes,
- shpFcn_);
- interscaleManager_.add_dense_matrix(nodalAtomicCount_,
- "AtomicDimensionlessMassMat");
- }
-#endif
+
for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
(_tiIt_->second)->construct_transfers();
}
@@ -293,24 +202,7 @@ namespace ATC {
// nothing needed in other cases since kinetostat force is balanced by boundary flux in FE equations
atomicRegulator_->reset_lambda_contribution(nodalAtomicFieldsRoc_[VELOCITY].quantity());
}
-#ifdef OBSOLETE
- //---------------------------------------------------------
- // compute_md_mass_matrix
- // compute the mass matrix arising from only atomistic
- // quadrature and contributions as a summation
- //---------------------------------------------------------
- void ATC_CouplingMomentum::compute_md_mass_matrix(FieldName thisField,
- DIAG_MAT & massMat)
- {
-
-
- if (thisField == DISPLACEMENT || thisField == VELOCITY)
- massMat.reset(nodalAtomicMass_->quantity());
- else if (thisField == MASS_DENSITY) { // dimensionless mass matrix
- massMat.reset(nodalAtomicCount_->quantity());
- }
- }
-#endif
+
//--------------------------------------------------------
// modify
// parses inputs and modifies state of the filter
diff --git a/lib/atc/ATC_CouplingMomentum.h b/lib/atc/ATC_CouplingMomentum.h
index d6ade453353493dd47116dd695ea814e9b6b74de..66a4b22ca75e979f9d13faf3c17b673ea58092be 100644
--- a/lib/atc/ATC_CouplingMomentum.h
+++ b/lib/atc/ATC_CouplingMomentum.h
@@ -44,14 +44,6 @@ namespace ATC {
/** pre time integration */
virtual void initialize();
- /** flags whether a methods reset is required */
- virtual bool reset_methods() const {
- bool resetMethods = ATC_Method::reset_methods() ||atomicRegulator_->need_reset();
- _ctiIt_ = timeIntegrators_.find(VELOCITY);
- if (_ctiIt_ == timeIntegrators_.end()) return resetMethods;
- return resetMethods || (_ctiIt_->second)->need_reset();
- };
-
/** pre/post atomic force calculation in minimize */
virtual void min_pre_force();
virtual void min_post_force();
@@ -76,17 +68,7 @@ namespace ATC {
//virtual void construct_time_integration_data();
/** set up data which is dependency managed */
virtual void construct_transfers();
-#ifdef OBSOLETE
- /** compute the mass matrix components coming from MD integration */
- virtual void compute_md_mass_matrix(FieldName thisField,
- DIAG_MAT & massMats);
- //
- /** operator to compute the mass matrix for the momentum equation from MD integration */
- AtfShapeFunctionRestriction * nodalAtomicMass_;
-
- /** operator to compute the dimensionless mass matrix from MD integration */
- AtfShapeFunctionRestriction * nodalAtomicCount_;
-#endif
+
/** physics specific filter initialization */
void init_filter();
diff --git a/lib/atc/ATC_CouplingMomentumEnergy.cpp b/lib/atc/ATC_CouplingMomentumEnergy.cpp
index 6ffaed5213dc2ec639af1b9ab40683ddd08c282b..7c1de1bbcb8a767881d3f09734e2adc3b69741c0 100644
--- a/lib/atc/ATC_CouplingMomentumEnergy.cpp
+++ b/lib/atc/ATC_CouplingMomentumEnergy.cpp
@@ -31,11 +31,6 @@ namespace ATC {
string matParamFile,
ExtrinsicModelType extrinsicModel)
: ATC_Coupling(groupName,perAtomArray,thisFix),
-#ifdef OBSOLETE
- nodalAtomicMass_(NULL),
- nodalAtomicCount_(NULL),
- nodalAtomicHeatCapacity_(NULL),
-#endif
nodalAtomicKineticTemperature_(NULL),
nodalAtomicConfigurationalTemperature_(NULL),
refPE_(0)
@@ -82,9 +77,6 @@ namespace ATC {
extVector_ = 1;
if (extrinsicModel != NO_MODEL)
sizeVector_ += extrinsicModelManager_.size_vector(sizeVector_);
-
- // create PE per atom ccompute
- lammpsInterface_->create_compute_pe_peratom();
}
//--------------------------------------------------------
@@ -111,88 +103,6 @@ namespace ATC {
// Base class initalizations
ATC_Coupling::initialize();
-
- // resetting precedence:
- // time integrator -> kinetostat/thermostat -> time filter
- // init_filter uses fieldRateNdFiltered which comes from the time integrator,
- // which is why the time integrator is initialized first
-
- // other initializations
-#ifdef OBSOLETE
-
- if (reset_methods()) {
- for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
- (_tiIt_->second)->force_reset();
- }
- atomicRegulator_->force_reset();
- }
-#endif
- if (reset_methods()) {
- for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
- (_tiIt_->second)->initialize();
- }
- atomicRegulator_->initialize();
- extrinsicModelManager_.initialize();
- }
- if (timeFilterManager_.need_reset()) {// reset thermostat power
- init_filter();
- }
- // clears need for reset
- timeFilterManager_.initialize();
- ghostManager_.initialize();
-
- if (!initialized_) {
- // initialize sources based on initial FE temperature
- double dt = lammpsInterface_->dt();
- prescribedDataMgr_->set_sources(time()+0.5*dt,sources_);
- extrinsicModelManager_.set_sources(fields_,extrinsicSources_);
- atomicRegulator_->compute_boundary_flux(fields_);
- compute_atomic_sources(fieldMask_,fields_,atomicSources_);
-
- // read in field data if necessary
- if (useRestart_) {
- RESTART_LIST data;
- read_restart_data(restartFileName_,data);
- useRestart_ = false;
- }
-
- // set consistent initial conditions, if requested
- if (!timeFilterManager_.filter_dynamics()) {
- if (consistentInitialization_) {
-
- DENS_MAT & velocity(fields_[VELOCITY].set_quantity());
- DENS_MAN * nodalAtomicVelocity(interscaleManager_.dense_matrix("NodalAtomicVelocity"));
- const DENS_MAT & atomicVelocity(nodalAtomicVelocity->quantity());
- DENS_MAT & temperature(fields_[TEMPERATURE].set_quantity());
- DENS_MAN * nodalAtomicTemperature(interscaleManager_.dense_matrix("NodalAtomicTemperature"));
- const DENS_MAT & atomicTemperature(nodalAtomicTemperature->quantity());
- const INT_ARRAY & nodeType(nodalGeometryType_->quantity());
- for (int i = 0; i<nNodes_; ++i) {
-
- if (nodeType(i,0)==MD_ONLY) {
- for (int j = 0; j < nsd_; j++) {
- velocity(i,j) = atomicVelocity(i,j);
- }
- temperature(i,0) = atomicTemperature(i,0);
- }
- }
- if (trackDisplacement_) {
- DENS_MAT & displacement(fields_[DISPLACEMENT].set_quantity());
- DENS_MAN * nodalAtomicDisplacement(interscaleManager_.dense_matrix("NodalAtomicDisplacement"));
- const DENS_MAT & atomicDisplacement(nodalAtomicDisplacement->quantity());
- for (int i = 0; i<nNodes_; ++i) {
-
- if (nodeType(i,0)==MD_ONLY) {
- for (int j = 0; j < nsd_; j++) {
- displacement(i,j) = atomicDisplacement(i,j);
- }
- }
- }
- }
- }
- }
- initialized_ = true;
- }
// reset integration field mask
intrinsicMask_.reset(NUM_FIELDS,NUM_FLUX);
@@ -368,38 +278,7 @@ namespace ATC {
TEMPERATURE);
interscaleManager_.add_dense_matrix(nodalAtomicTemperature,
"NodalAtomicTemperature");
-#ifdef OBSOLETE
- if (!useFeMdMassMatrix_) {
- // atomic momentum mass matrix
- FundamentalAtomQuantity * atomicMass = interscaleManager_.fundamental_atom_quantity(LammpsInterface::ATOM_MASS);
- nodalAtomicMass_ = new AtfShapeFunctionRestriction(this,
- atomicMass,
- shpFcn_);
- interscaleManager_.add_dense_matrix(nodalAtomicMass_,
- "AtomicMomentumMassMat");
-
- // atomic dimensionless mass matrix
- ConstantQuantity<double> * atomicOnes = new ConstantQuantity<double>(this,1);
- interscaleManager_.add_per_atom_quantity(atomicOnes,"AtomicOnes");
- nodalAtomicCount_ = new AtfShapeFunctionRestriction(this,
- atomicOnes,
- shpFcn_);
- interscaleManager_.add_dense_matrix(nodalAtomicCount_,
- "AtomicDimensionlessMassMat");
-
- // classical thermodynamic heat capacity of the atoms
- HeatCapacity * heatCapacity = new HeatCapacity(this);
- interscaleManager_.add_per_atom_quantity(heatCapacity,
- "AtomicHeatCapacity");
-
- // atomic thermal mass matrix
- nodalAtomicHeatCapacity_ = new AtfShapeFunctionRestriction(this,
- heatCapacity,
- shpFcn_);
- interscaleManager_.add_dense_matrix(nodalAtomicHeatCapacity_,
- "NodalAtomicHeatCapacity");
- }
-#endif
+
for (_tiIt_ = timeIntegrators_.begin(); _tiIt_ != timeIntegrators_.end(); ++_tiIt_) {
(_tiIt_->second)->construct_transfers();
}
@@ -431,28 +310,7 @@ namespace ATC {
atomicRegulator_->reset_lambda_contribution(powerMat,TEMPERATURE);
}
}
-#ifdef OBSOLETE
- //---------------------------------------------------------
- // compute_md_mass_matrix
- // compute the mass matrix arising from only atomistic
- // quadrature and contributions as a summation
- //---------------------------------------------------------
- void ATC_CouplingMomentumEnergy::compute_md_mass_matrix(FieldName thisField,
- DIAG_MAT & massMat)
- {
-
-
- if (thisField == DISPLACEMENT || thisField == VELOCITY) {
- massMat.reset(nodalAtomicMass_->quantity());
- }
- else if (thisField == MASS_DENSITY) { // dimensionless mass matrix
- massMat.reset(nodalAtomicCount_->quantity());
- }
- else if (thisField == TEMPERATURE) {
- massMat.reset(nodalAtomicHeatCapacity_->quantity());
- }
- }
-#endif
+
//--------------------------------------------------------
// modify
// parses inputs and modifies state of the filter
diff --git a/lib/atc/ATC_CouplingMomentumEnergy.h b/lib/atc/ATC_CouplingMomentumEnergy.h
index d0c78bdaeea971574cc1e8b8e7d6eb23910c5cdb..598ee13bd7bea83fd28d8d3083bdddc683dfeb8b 100644
--- a/lib/atc/ATC_CouplingMomentumEnergy.h
+++ b/lib/atc/ATC_CouplingMomentumEnergy.h
@@ -42,16 +42,6 @@ namespace ATC {
/** pre time integration */
virtual void initialize();
- /** flags whether a methods reset is required */
-
- virtual bool reset_methods() const {
- bool resetMethods = ATC_Method::reset_methods() || atomicRegulator_->need_reset();
- for (_ctiIt_ = timeIntegrators_.begin(); _ctiIt_ != timeIntegrators_.end(); ++_ctiIt_) {
- resetMethods |= (_ctiIt_->second)->need_reset();
- }
- return resetMethods;
- };
-
/** compute scalar for output - added energy */
virtual double compute_scalar(void);
@@ -72,20 +62,7 @@ namespace ATC {
//virtual void construct_time_integration_data();
/** set up data which is dependency managed */
virtual void construct_transfers();
-#ifdef OBSOLETE
- /** compute the mass matrix components coming from MD integration */
- virtual void compute_md_mass_matrix(FieldName thisField,
- DIAG_MAT & massMats);
- //
- /** operator to compute the mass matrix for the momentum equation from MD integration */
- AtfShapeFunctionRestriction * nodalAtomicMass_;
-
- /** operator to compute the dimensionless mass matrix from MD integration */
- AtfShapeFunctionRestriction * nodalAtomicCount_;
-
- /** operator to compute mass matrix from MD */
- AtfShapeFunctionRestriction * nodalAtomicHeatCapacity_;
-#endif
+
/** physics specific filter initialization */
void init_filter();
/** kinetic temperature for post-processing */
diff --git a/lib/atc/ATC_Method.cpp b/lib/atc/ATC_Method.cpp
index 90d29f853edf923247588f7edd14faa78ed45ca9..f32b06cb6ff35d796dc1a3bb85a298a56f0fd6aa 100644
--- a/lib/atc/ATC_Method.cpp
+++ b/lib/atc/ATC_Method.cpp
@@ -699,6 +699,8 @@ pecified
<TT> fix_modify AtC internal_atom_integrate on </TT>
\section description
Has AtC perform time integration for the atoms in the group on which it operates. This does not include boundary atoms.
+ \section restrictions
+ AtC must be created before any fixes doing time integration. It must be on for coupling methods which impose constraints on velocities during the first verlet step, e.g. control momentum glc_velocity.
\section default
on for coupling methods, off for post-processors
off
@@ -1503,16 +1505,17 @@ pecified
}
}
//--------------------------------------------------------
- void ATC_Method::init_integrate_velocity()
+ void ATC_Method::init_integrate()
{
atomTimeIntegrator_->init_integrate_velocity(dt());
ghostManager_.init_integrate_velocity(dt());
- }
- //--------------------------------------------------------
- void ATC_Method::init_integrate_position()
- {
+ // account for other fixes doing time integration
+ interscaleManager_.fundamental_force_reset(LammpsInterface::ATOM_VELOCITY);
+
atomTimeIntegrator_->init_integrate_position(dt());
ghostManager_.init_integrate_position(dt());
+ // account for other fixes doing time integration
+ interscaleManager_.fundamental_force_reset(LammpsInterface::ATOM_POSITION);
}
//-------------------------------------------------------------------
void ATC_Method::post_init_integrate()
@@ -1564,6 +1567,8 @@ pecified
{
atomTimeIntegrator_->final_integrate(dt());
ghostManager_.final_integrate(dt());
+ // account for other fixes doing time integration
+ interscaleManager_.fundamental_force_reset(LammpsInterface::ATOM_VELOCITY);
}
//-------------------------------------------------------------------
void ATC_Method::post_final_integrate()
@@ -1577,7 +1582,6 @@ pecified
{
localStep_ += 1;
}
-
//--------------------------------------------------------------
void ATC_Method::finish()
{
diff --git a/lib/atc/ATC_Method.h b/lib/atc/ATC_Method.h
index 17cd859dc1276990def9f9809c5d080f94c2e251..d7708d298874b6e04af07dca3707b69d34fc5e62 100644
--- a/lib/atc/ATC_Method.h
+++ b/lib/atc/ATC_Method.h
@@ -68,14 +68,8 @@ namespace ATC {
update_step();
};
- /** Predictor phase, Verlet first step for velocity */
- virtual void init_integrate_velocity();
-
- /** Predictor phase, executed between velocity and position Verlet */
- virtual void mid_init_integrate(){};
-
- /** Predictor phase, Verlet first step for position */
- virtual void init_integrate_position();
+ /** Predictor phase, Verlet first step for velocity and position */
+ virtual void init_integrate();
/** Predictor phase, executed after Verlet */
virtual void post_init_integrate();
diff --git a/lib/atc/ATC_Transfer.cpp b/lib/atc/ATC_Transfer.cpp
index 0e47d5fc1b9a3aacd5e6ec4198b38330698d26f9..968dcbbea71ae3c610d0ff486e424617a8482a12 100644
--- a/lib/atc/ATC_Transfer.cpp
+++ b/lib/atc/ATC_Transfer.cpp
@@ -911,6 +911,8 @@ namespace ATC {
if ( output_now() && !outputStepZero_ ) output();
outputStepZero_ = false;
+ //ATC_Method::post_final_integrate();
+
}
//-------------------------------------------------------------------
diff --git a/lib/atc/ATC_Transfer.h b/lib/atc/ATC_Transfer.h
index 8f4753bab72df9bfcc4ef37948f00fcfe7e267dd..9e05b8b84bcaac7401d9b2a50dafc31a3a2611c3 100644
--- a/lib/atc/ATC_Transfer.h
+++ b/lib/atc/ATC_Transfer.h
@@ -46,7 +46,6 @@ class ATC_Transfer : public ATC_Method {
/** second time substep routine */
virtual void pre_final_integrate();
- //virtual void final_integrate(){};
virtual void post_final_integrate();
/** communication routines */
diff --git a/lib/atc/CbEam.h b/lib/atc/CbEam.h
index 1ff7178c004941efb40dceb2be7e34eaa2397a57..7b7df175e825094900f59c76eeb77d19378b83dd 100644
--- a/lib/atc/CbEam.h
+++ b/lib/atc/CbEam.h
@@ -9,7 +9,6 @@
namespace ATC
{
- class LAMMPS_NS::PairEAM; // necessary for non-lib build of ATC (vs USER-ATC)
/**
* @class CbEam
diff --git a/lib/atc/ConcentrationRegulator.cpp b/lib/atc/ConcentrationRegulator.cpp
index 430c3a7d41fd18e581f8ca76eade56fb5ec96562..c734a5d98db3b49aea427223d398fc54708e6579 100644
--- a/lib/atc/ConcentrationRegulator.cpp
+++ b/lib/atc/ConcentrationRegulator.cpp
@@ -562,7 +562,7 @@ Tv(0) = 300.;
int ConcentrationRegulatorMethodTransition::pick_element() const
{
double r = uniform();
- ESET::iterator itr = elemset_.begin(); // global?
+ ESET::const_iterator itr = elemset_.begin(); // global?
for (int i = 0; i < volumes_.size() ; ++i) {
if (r < volumes_(i)) return *itr;
itr++;
diff --git a/lib/atc/ElasticTimeIntegrator.cpp b/lib/atc/ElasticTimeIntegrator.cpp
index b8058192efec719872b7335064ee46e09c69cda6..c2a592d46c73af0907ac531f82417cf1e21da32a 100644
--- a/lib/atc/ElasticTimeIntegrator.cpp
+++ b/lib/atc/ElasticTimeIntegrator.cpp
@@ -230,10 +230,10 @@ namespace ATC {
}
//--------------------------------------------------------
- // mid_initial_integrate1
- // time integration at the mid-point of Verlet step 1
+ // pre_initial_integrate1
+ // time integration before Verlet step 1
//--------------------------------------------------------
- void ElasticTimeIntegratorVerlet::mid_initial_integrate1(double dt)
+ void ElasticTimeIntegratorVerlet::pre_initial_integrate1(double dt)
{
explicit_1(velocity_.set_quantity(),acceleration_.quantity(),.5*dt);
}
@@ -334,10 +334,10 @@ namespace ATC {
}
//--------------------------------------------------------
- // mid_initial_integrate1
- // time integration at the mid-point of Verlet step 1
+ // pre_initial_integrate1
+ // time integration before Verlet step 1
//--------------------------------------------------------
- void ElasticTimeIntegratorVerletFiltered::mid_initial_integrate1(double dt)
+ void ElasticTimeIntegratorVerletFiltered::pre_initial_integrate1(double dt)
{
explicit_1(velocity_.set_quantity(),acceleration_.quantity(),.5*dt);
explicit_1(nodalAtomicVelocityOut_.set_quantity(),nodalAtomicAcceleration_.quantity(),.5*dt);
@@ -497,11 +497,10 @@ namespace ATC {
}
//--------------------------------------------------------
- // mid_initial_integrate1
- // time integration between the velocity and position
- // updates in Verlet step 1
+ // post_initial_integrate1
+ // time integration after Verlet step 1
//--------------------------------------------------------
- void ElasticTimeIntegratorFractionalStep::mid_initial_integrate1(double dt)
+ void ElasticTimeIntegratorFractionalStep::post_initial_integrate1(double dt)
{
// atomic contributions to change in momentum
// compute change in restricted atomic momentum
@@ -531,16 +530,9 @@ namespace ATC {
nodalAtomicVelocityOld_,
VELOCITY);
nodalAtomicMomentumOld_ = nodalAtomicMomentum;
- }
- //--------------------------------------------------------
- // post_initial_integrate1
- // time integration after Verlet step 1
- //--------------------------------------------------------
- void ElasticTimeIntegratorFractionalStep::post_initial_integrate1(double dt)
- {
// get nodal momentum for second part of force update
- nodalAtomicForce_ = nodalAtomicMomentum_->quantity();
+ nodalAtomicForce_ = nodalAtomicMomentum;
nodalAtomicForce_ *= -1.;
// update nodal displacements
diff --git a/lib/atc/ElasticTimeIntegrator.h b/lib/atc/ElasticTimeIntegrator.h
index f6d270f5fa19bc3f2e9852d56a3d1fb2df013b10..02e6b5c3c860377d001dda1e268e95025cf6e415 100644
--- a/lib/atc/ElasticTimeIntegrator.h
+++ b/lib/atc/ElasticTimeIntegrator.h
@@ -131,8 +131,8 @@ namespace ATC {
virtual void initialize();
// time step methods, corresponding to ATC_Transfer
- /** first part of mid_initial_integrate */
- virtual void mid_initial_integrate1(double dt);
+ /** first part of pre_initial_integrate */
+ virtual void pre_initial_integrate1(double dt);
/** first part of post_initial_integrate */
virtual void post_initial_integrate1(double dt);
/** first part of pre_final_integrate */
@@ -190,8 +190,8 @@ namespace ATC {
virtual ~ElasticTimeIntegratorVerletFiltered(){};
// time step methods, corresponding to ATC_Transfer
- /** first part of mid_initial_integrate */
- virtual void mid_initial_integrate1(double dt);
+ /** first part of pre_initial_integrate */
+ virtual void pre_initial_integrate1(double dt);
/** first part of post_initial_integrate */
virtual void post_initial_integrate1(double dt);
/** second part of post_final_integrate */
@@ -243,10 +243,9 @@ namespace ATC {
// time step methods, corresponding to ATC_Transfer
/** first part of pre_initial_integrate */
virtual void pre_initial_integrate1(double dt);
- /** second part of mid_initial_integrate */
+ /** second part of pre_initial_integrate */
virtual void pre_initial_integrate2(double dt);
- /** first part of mid_initial_integrate */
- virtual void mid_initial_integrate1(double dt);
+
/** first part of post_initial_integrate */
virtual void post_initial_integrate1(double dt);
/** second part of post_final_integrate */
diff --git a/lib/atc/ExtrinsicModel.cpp b/lib/atc/ExtrinsicModel.cpp
index 0522f068193f6a7f198e331896c0d931019fce48..cdf556f4c1a1f5b7f98bf859efbf4ad5bad39285 100644
--- a/lib/atc/ExtrinsicModel.cpp
+++ b/lib/atc/ExtrinsicModel.cpp
@@ -233,23 +233,6 @@ namespace ATC {
}
}
- //--------------------------------------------------------
- // mid_init_integrate
- //--------------------------------------------------------
- void ExtrinsicModelManager::mid_init_integrate(ExtrinsicModelType modelType)
- {
- vector<ExtrinsicModel *>::iterator imodel;
- if (modelType == NUM_MODELS) {// execute all the models
- for(imodel=extrinsicModels_.begin(); imodel!=extrinsicModels_.end(); imodel++)
- (*imodel)->mid_init_integrate();
- }
- else { // execute only requested type of model
- for(imodel=extrinsicModels_.begin(); imodel!=extrinsicModels_.end(); imodel++)
- if ((*imodel)->model_type() == modelType)
- (*imodel)->mid_init_integrate();
- }
- }
-
//--------------------------------------------------------
// post_init_integrate
//--------------------------------------------------------
diff --git a/lib/atc/ExtrinsicModel.h b/lib/atc/ExtrinsicModel.h
index 4d8fb95fd782485685476777790d5ff52e6e7de7..43ad08d21139c848e7d0d111babd3b8e465eb3fc 100644
--- a/lib/atc/ExtrinsicModel.h
+++ b/lib/atc/ExtrinsicModel.h
@@ -77,8 +77,7 @@ namespace ATC {
// calls during LAMMPS Velocity-Verlet integration
/** Predictor phase, executed before Verlet */
void pre_init_integrate(ExtrinsicModelType modelType = NUM_MODELS);
- /** Predictor phase, executed between velocity and position Verlet */
- void mid_init_integrate(ExtrinsicModelType modelType = NUM_MODELS);
+
/** Predictor phase, executed after Verlet */
void post_init_integrate(ExtrinsicModelType modelType = NUM_MODELS);
@@ -244,9 +243,6 @@ namespace ATC {
/** Predictor phase, executed before Verlet */
virtual void pre_init_integrate(){};
- /** Predictor phase, executed between velocity and position Verlet */
- virtual void mid_init_integrate(){};
-
/** Predictor phase, executed after Verlet */
virtual void post_init_integrate(){};
diff --git a/lib/atc/FE_Engine.cpp b/lib/atc/FE_Engine.cpp
index 14bcb58c8b38e2276ee59cfd6651378648d0e5a5..a87a03032cda4d9369c0be48099fc6282b60f3d7 100644
--- a/lib/atc/FE_Engine.cpp
+++ b/lib/atc/FE_Engine.cpp
@@ -89,7 +89,7 @@ namespace ATC{
// now do all FE_Engine data structure partitioning
// partition nullElements_
- /*for (vector<int>::iterator elemsIter = feMesh_->processor_elts().begin();
+ /*for (vector<int>::const_iterator elemsIter = feMesh_->processor_elts().begin();
elemsIter != feMesh_->processor_elts().end();
++elemsIter)
{
@@ -550,7 +550,7 @@ namespace ATC{
DENS_MAT elementMassMatrix(nNodesPerElement_,nNodesPerElement_);
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -604,7 +604,7 @@ namespace ATC{
DENS_MAT coefsAtIPs;
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -782,7 +782,7 @@ namespace ATC{
// element wise assembly
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -849,7 +849,7 @@ namespace ATC{
massMatrix.reset(nNodesUnique_,nNodesUnique_);// zero since partial fill
DENS_MAT elementMassMatrix(nNodesPerElement_,nNodesPerElement_);
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -905,7 +905,7 @@ namespace ATC{
// assemble lumped diagonal mass
DENS_VEC Nvec(nNodesPerElement_);
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -944,7 +944,7 @@ namespace ATC{
}
// assemble diagonal matrix
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -1143,7 +1143,7 @@ namespace ATC{
Array<int> count(nNodesUnique_); count = 0;
set_quadrature(NODAL);
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -1204,7 +1204,7 @@ namespace ATC{
DENS_MAT elementEnergy(nNodesPerElement_,1); // workspace
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -1269,7 +1269,7 @@ namespace ATC{
// Iterate over elements partitioned onto current processor.
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -1337,7 +1337,7 @@ namespace ATC{
}
}
- vector<FieldName>::iterator _fieldIter_;
+ vector<FieldName>::const_iterator _fieldIter_;
for (_fieldIter_ = usedFields.begin(); _fieldIter_ != usedFields.end();
++_fieldIter_) {
// myRhs is where we put the global result for this field.
@@ -1660,7 +1660,7 @@ namespace ATC{
// element wise assembly
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -1740,7 +1740,7 @@ namespace ATC{
}
// element wise assembly
- set< PAIR >::iterator iter;
+ set< PAIR >::const_iterator iter;
for (iter = faceSet.begin(); iter != faceSet.end(); iter++) {
// get connectivity
int ielem = iter->first;
@@ -2214,7 +2214,7 @@ namespace ATC{
}
vector<int> myElems = feMesh_->owned_elts();
- for (vector<int>::iterator elemsIter = myElems.begin();
+ for (vector<int>::const_iterator elemsIter = myElems.begin();
elemsIter != myElems.end();
++elemsIter)
{
@@ -2291,7 +2291,7 @@ namespace ATC{
// the data member?
// element wise assembly
- set< pair <int,int> >::iterator iter;
+ set< pair <int,int> >::const_iterator iter;
for (iter = faceSet.begin(); iter != faceSet.end(); iter++) {
int ielem = iter->first;
// if this is not our element, do not do calculations
diff --git a/lib/atc/FE_Mesh.cpp b/lib/atc/FE_Mesh.cpp
index 3dd00497d6f9e5c03695bf2c6b8f0b1a9eb4bbdd..e97755bcf19e4ebd73dc8dc18dc430d634c0c2be 100644
--- a/lib/atc/FE_Mesh.cpp
+++ b/lib/atc/FE_Mesh.cpp
@@ -1432,7 +1432,8 @@ namespace ATC {
int FE_Mesh::map_elem_to_myElem(int elemID) const
{
- return elemToMyElemMap_.at(elemID);
+
+ return elemToMyElemMap_.find(elemID)->second;
}
int FE_Mesh::map_myElem_to_elem(int myElemID) const
diff --git a/lib/atc/InterscaleOperators.cpp b/lib/atc/InterscaleOperators.cpp
index 28762b9772f648a4b995f641e9db48c23d9a4f2d..0e218bdef438aed3190dc01260779fdc9cf2d4cd 100644
--- a/lib/atc/InterscaleOperators.cpp
+++ b/lib/atc/InterscaleOperators.cpp
@@ -507,15 +507,23 @@ namespace ATC{
//--------------------------------------------------------
// computes_force_reset
//--------------------------------------------------------
- void InterscaleManager::lammps_force_reset()
+ void InterscaleManager::fundamental_force_reset(unsigned quantity)
{
for (unsigned int i = 0; i < NUM_ATOM_TYPES; i++) {
- for (unsigned int j = 0; j < LammpsInterface::NUM_FUNDAMENTAL_ATOM_QUANTITIES; j++) {
- if (fundamentalAtomQuantities_[i][j]) {
- fundamentalAtomQuantities_[i][j]->lammps_force_reset();
- }
+ if (fundamentalAtomQuantities_[i][quantity]) {
+ fundamentalAtomQuantities_[i][quantity]->lammps_force_reset();
}
}
+ }
+
+ //--------------------------------------------------------
+ // computes_force_reset
+ //--------------------------------------------------------
+ void InterscaleManager::lammps_force_reset()
+ {
+ for (unsigned int j = 0; j < LammpsInterface::NUM_FUNDAMENTAL_ATOM_QUANTITIES; j++) {
+ fundamental_force_reset(j);
+ }
lammps_reset_loop(perAtomQuantities_);
lammps_reset_loop(perAtomIntQuantities_);
lammps_reset_loop(perAtomDiagonalMatrices_);
diff --git a/lib/atc/InterscaleOperators.h b/lib/atc/InterscaleOperators.h
index 20d16de9c502fe91f0bf8e9f693d888636cb76a5..4681801e732dfacf1e2510cff8e09d70b08d9bd5 100644
--- a/lib/atc/InterscaleOperators.h
+++ b/lib/atc/InterscaleOperators.h
@@ -161,6 +161,9 @@ namespace ATC {
/** resets nlocal count of managed atomic quantities which do not perform parallel exchange */
void reset_nlocal();
+ /** resets specific lammps fundamental quantities data, as needed, to account for times when lammps can change quantities */
+ void fundamental_force_reset(unsigned quantity);
+
/** resets all lammps data, as needed, to account for times when lammps can change quantities */
void lammps_force_reset();
diff --git a/lib/atc/Kinetostat.cpp b/lib/atc/Kinetostat.cpp
index a76c678d220ac726e4162dc66a603995fe987e3f..c104efc720f183e8a6659299680723f8b908c0d7 100644
--- a/lib/atc/Kinetostat.cpp
+++ b/lib/atc/Kinetostat.cpp
@@ -1129,7 +1129,7 @@ namespace ATC {
// apply the kinetostat to the atoms in the
// mid-predictor integration phase
//--------------------------------------------------------
- void StressFlux::apply_mid_predictor(double dt)
+ void StressFlux::apply_pre_predictor(double dt)
{
double dtLambda = 0.5*dt;
// apply lambda force to atoms
diff --git a/lib/atc/Kinetostat.h b/lib/atc/Kinetostat.h
index e85caa8860466121d1ffabdb82c1e757ea0f11a0..18400342962dfe67ce3937ff936462badd5c66d4 100644
--- a/lib/atc/Kinetostat.h
+++ b/lib/atc/Kinetostat.h
@@ -389,8 +389,8 @@ namespace ATC {
/** instantiate all needed data */
virtual void construct_transfers();
- /** applies kinetostat to atoms in the mid-predictor phase */
- virtual void apply_mid_predictor(double dt);
+ /** applies kinetostat to atoms in the pre-predictor phase */
+ virtual void apply_pre_predictor(double dt);
/** applies kinetostat to atoms in the post-corrector phase */
virtual void apply_post_corrector(double dt);
diff --git a/lib/atc/SpeciesTimeIntegrator.h b/lib/atc/SpeciesTimeIntegrator.h
index 5d6c962d719336373935d32f9d193f1206783d06..c8fbe97c5e83b90b72267e0ee3f1efe2418b0a36 100644
--- a/lib/atc/SpeciesTimeIntegrator.h
+++ b/lib/atc/SpeciesTimeIntegrator.h
@@ -84,7 +84,7 @@ namespace ATC {
DENS_MAN & massDensity_;
/** atomic nodal mass density field */
-// OBOLETE?
+// OBSOLETE?
DENS_MAN & nodalAtomicMassDensityOut_;
DENS_MAN * nodalAtomicMassDensity_;
diff --git a/lib/atc/TimeIntegrator.cpp b/lib/atc/TimeIntegrator.cpp
index 26db421b5ca434f693669fae80fc60af06afc10d..4e8dec8c6299f256e0a85a6f121aefc4e1bbb1db 100644
--- a/lib/atc/TimeIntegrator.cpp
+++ b/lib/atc/TimeIntegrator.cpp
@@ -45,7 +45,7 @@ namespace ATC {
void AtomTimeIntegratorType::init_integrate_velocity(double dt)
{
const DENS_MAT & m(mass_->quantity());
-
+
_deltaQuantity_ = force_->quantity();
_deltaQuantity_ /= m;
_deltaQuantity_ *= 0.5*dt;
@@ -148,26 +148,6 @@ namespace ATC {
timeIntegrationMethod_->pre_initial_integrate2(dt);
}
- //--------------------------------------------------------
- // mid_initial_integrate1
- // first time integration computations
- // at the mid-point of Verlet step 1
- //--------------------------------------------------------
- void TimeIntegrator::mid_initial_integrate1(double dt)
- {
- timeIntegrationMethod_->mid_initial_integrate1(dt);
- }
-
- //--------------------------------------------------------
- // mid_initial_integrate2
- // second time integration computations
- // at the mid-point of Verlet step 1
- //--------------------------------------------------------
- void TimeIntegrator::mid_initial_integrate2(double dt)
- {
- timeIntegrationMethod_->mid_initial_integrate2(dt);
- }
-
//--------------------------------------------------------
// post_initial_integrate1
// first time integration computations
diff --git a/lib/atc/TimeIntegrator.h b/lib/atc/TimeIntegrator.h
index 5c29a6e121a2b50db20e4fffe4c354f2d8713be7..22e4b716fd771048f76678b57e41b7543598144e 100644
--- a/lib/atc/TimeIntegrator.h
+++ b/lib/atc/TimeIntegrator.h
@@ -145,12 +145,7 @@ namespace ATC {
virtual void pre_initial_integrate1(double dt);
/** second part of pre_initial_integrate */
virtual void pre_initial_integrate2(double dt);
-
- /** first part of mid_initial_integrate */
- virtual void mid_initial_integrate1(double dt);
- /** second part of mid_initial_integrate */
- virtual void mid_initial_integrate2(double dt);
-
+
/** first part of post_initial_integrate */
virtual void post_initial_integrate1(double dt);
/** second part of post_initial_integrate */
@@ -257,12 +252,7 @@ namespace ATC {
virtual void pre_initial_integrate1(double dt){};
/** second part of pre_initial_integrate */
virtual void pre_initial_integrate2(double dt){};
-
- /** first part of mid_initial_integrate */
- virtual void mid_initial_integrate1(double dt){};
- /** second part of mid_initial_integrate */
- virtual void mid_initial_integrate2(double dt){};
-
+
/** first part of post_initial_integrate */
virtual void post_initial_integrate1(double dt){};
/** second part of post_initial_integrate */