diff --git a/src/common/action_interface.cc b/src/common/action_interface.cc index 5142bbc..bbf7d7f 100644 --- a/src/common/action_interface.cc +++ b/src/common/action_interface.cc @@ -1,202 +1,202 @@ /** * @file action_interface.cc * * @author Guillaume Anciaux * @author Till Junge * * @date Mon Sep 08 23:40:22 2014 * * @brief This is the interface to all action objects * * @section LICENSE * * Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * LibMultiScale is free software: you can redistribute it and/or modify it * under the * terms of the GNU Lesser General Public License as published by the Free * Software Foundation, either version 3 of the License, or (at your option) any * later version. * * LibMultiScale is distributed in the hope that it will be useful, but * WITHOUT ANY * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR * A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more * details. * * You should have received a copy of the GNU Lesser General Public License * along with LibMultiScale. If not, see . * */ #include "action_interface.hh" #include "filter_interface.hh" #include "lm_common.hh" __BEGIN_LIBMULTISCALE__ /* -------------------------------------------------------------------------- */ ActionInterface::ActionInterface() { setDefaults(); } /* -------------------------------------------------------------------------- */ void ActionInterface::printself(std::ostream &stream) const { stream << "Action " << this->getID() << " " << typeid(*this).name() << "action_step " << action_step << " " << "start_step " << start_step << " " << "end_step " << end_step << " " << "frequency " << frequency << std::endl; } /* -------------------------------------------------------------------------- */ void ActionInterface::setDefaults() { action_step = 0; frequency = 1; start_step = 0; end_step = UINT_MAX; one_shot = UINT_MAX; } /* -------------------------------------------------------------------------- */ bool ActionInterface::shouldMakeAction() { if ((current_stage == PRE_FATAL) && (stage_mask & PRE_FATAL)) return true; bool sM = doesStageMaskMatch(); if (!sM) return false; setOneStep(); DUMP("for action " << this->getID() << " freq " << frequency << " step " << action_step, DBG_INFO); if (end_step >= current_step && start_step <= current_step) if (frequency != UINT_MAX && current_step % frequency == 0) { DUMP("should proceed action " << this->getID() << " current_step = " << current_step << " start " << start_step << " end " << end_step, DBG_INFO); return true; } return false; } /* -------------------------------------------------------------------------- */ void ActionInterface::setOneStep() { if (one_shot != UINT_MAX) end_step = start_step = one_shot; } /* -------------------------------------------------------------------------- */ bool ActionInterface::doesStageMaskMatch() { if (stage_mask.isNone()) stage_mask = PRE_STEP1; if (!(stage_mask & current_stage)) return false; return true; } /* -------------------------------------------------------------------------- */ void ActionInterface::action() { - this->compute(); + this->compute(true); ++action_step; } /* -------------------------------------------------------------------------- */ /* LMDESC ActionInterface This describes the generic keywords which all actions inheritate from. */ void ActionInterface::declareParams() { /* LMKEYWORD INPUT Specify the input(s) of the Action component */ this->parseKeyword("INPUT", static_cast(*this)); this->makeItOptional("INPUT"); /* LMKEYWORD FREQ Fix the frequency at which action should be called */ this->parseKeyword("FREQ", frequency, 1u); /* LMKEYWORD START Fix the first step at which action should start */ this->parseKeyword("START", start_step, 0u); /* LMKEYWORD END Fix the last step at which action should be performed */ this->parseKeyword("END", end_step, lm_uint_max); /* LMKEYWORD ONESHOT Set the action to be called only at a given step. */ this->parseKeyword("ONESHOT", one_shot, lm_uint_max); /* LMKEYWORD STAGE Set the stage within generic explicit integration scheme at which the action should be called. The possible values are: - PRE_DUMP - PRE_STEP1 - PRE_STEP2 - PRE_STEP3 - PRE_STEP4 - PRE_STEP5 - PRE_STEP6 - PRE_STEP7 The keyword STAGE can be called more than once and the result will be a bit-mask value so that an action can be called at several stages within a given timestep. These stages are identified depending on the implemented integration scheme. For instance the classical Verlet integration scheme is coded in a loop in the AMEL executable as follows: .. code-block:: cpp for (int current_step = 0; current_step < nb_step ; ++current_step, current_time += dt){ stimulator.stimulate(PRE_DUMP); dumper.dump(); stimulator.stimulate(PRE_STEP1); dom.performStep1(); dom.coupling(COUPLING_STEP1); stimulator.stimulate(PRE_STEP2); dom.performStep2(); dom.coupling(COUPLING_STEP2); stimulator.stimulate(PRE_STEP3); dom.performStep3(); dom.coupling(COUPLING_STEP3); stimulator.stimulate(PRE_STEP4); dom.coupling(COUPLING_STEP4); } For a classical Verlet integration scheme performStep1() is the initial velocity and position update, performStep2() is the force computation, and performStep3() is the velocity correction. If the stage value is kept to NONE, the default of PRE\_STEP1 will apply automatically */ this->parseKeyword("STAGE", stage_mask, NONE_STEP); } /* -------------------------------------------------------------------------- */ __END_LIBMULTISCALE__ diff --git a/src/common/component_libmultiscale.cc b/src/common/component_libmultiscale.cc index 6ce158f..17d7fb7 100644 --- a/src/common/component_libmultiscale.cc +++ b/src/common/component_libmultiscale.cc @@ -1,355 +1,369 @@ /* -------------------------------------------------------------------------- */ #include "component_libmultiscale.hh" #include "comm_group.hh" #include "container_array.hh" #include "lm_communicator.hh" /* -------------------------------------------------------------------------- */ __BEGIN_LIBMULTISCALE__ /* -------------------------------------------------------------------------- */ #define secure(...) secure_call([&]() { __VA_ARGS__; }); /* --------------------------------------------------------------------- */ // OutputContainer /* --------------------------------------------------------------------- */ std::ostream &operator<<(std::ostream &os, OutputContainer &out) { out.get().printself(os); return os; } OutputContainer &OutputContainer::operator=(Component &v) { AutoDispatch::ArgumentAny::operator=(v); this->component = &v; return *this; } UInt OutputContainer::evalRelease() const { UInt release = this->get().getRelease(); if (component) return hash_releases(release, component->evalRelease()); return release; } /* ---------------------------------------------------------------------- */ UInt OutputContainer::getRelease() const { UInt release = this->get().getRelease(); if (component) return hash_releases(release, component->getRelease()); return release; } /* --------------------------------------------------------------------- */ // InputContainer /* --------------------------------------------------------------------- */ InputContainer &InputContainer::operator=(OutputContainer &v) { output_reference = &v; return *this; } /* --------------------------------------------------------------------- */ InputContainer &InputContainer::operator=(OutputContainer &&v) { output_reference = std::make_shared(std::forward(v)); return *this; } OutputContainer &InputContainer::get() const { if (not output_reference.has_value()) throw Component::UnconnectedInput{"input is not connected"}; try { auto *ptr = std::any_cast(output_reference); return *ptr; } catch (std::bad_any_cast &e) { } auto ptr = std::any_cast>(output_reference); return *ptr; } /* --------------------------------------------------------------------- */ OutputContainer &InputContainer::eval() const { OutputContainer &out = this->get(); if (out.component != nullptr) { out.component->compute(); } return out; } /* --------------------------------------------------------------------- */ UInt InputContainer::evalRelease() const { OutputContainer &out = this->get(); return out.evalRelease(); } /* --------------------------------------------------------------------- */ UInt InputContainer::getRelease() const { OutputContainer &out = this->get(); return out.getRelease(); } /* -------------------------------------------------------------------------- */ // Component /* -------------------------------------------------------------------------- */ std::ostream &operator<<(std::ostream &os, Component &comp) { comp.printself(os); return os; } /* --------------------------------------------------------------------- */ Component::Component() : calculated_once(false) {} /* --------------------------------------------------------------------- */ InputContainer & Component::getInternalInput(const std::string &requested_input) { auto &&[input, arg_cont_input] = getContainer(requested_input, inputs); if (not arg_cont_input.has_value()) throw UnconnectedInput{"For component '" + this->getID() + "' input '" + input + "' is not connected"}; return arg_cont_input; } /* --------------------------------------------------------------------- */ OutputContainer &Component::getInput(const std::string &requested_input) { return this->getInternalInput(requested_input).eval(); } /* --------------------------------------------------------------------- */ OutputContainer &Component::getOutput(const std::string &output_name) { auto &&[name, output] = getContainer(output_name, outputs); return output; } /* --------------------------------------------------------------------- */ OutputContainer &Component::evalOutput(const std::string &output_name) { this->compute(); return getOutput(output_name); } /* --------------------------------------------------------------------- */ auto Component::evalOutputs() -> const decltype(outputs) & { this->compute(); return outputs; } /* --------------------------------------------------------------------- */ auto Component::getOutputs() -> const decltype(outputs) & { return outputs; } /* --------------------------------------------------------------------- */ void Component::printself(std::ostream &os) const { for (auto &pair : this->inputs) { os << "input " << pair.first << ":"; if (pair.second.has_value()) os << pair.second.get() << "\n"; else os << "not yet defined/computed\n"; } for (auto &pair : this->outputs) { os << "output " << pair.first << ":"; try { os << pair.second.get() << "\n"; } catch (...) { os << "not yet defined/computed\n"; } } } /* --------------------------------------------------------------------- */ void Component::connect(const std::string &input, Component &comp, const std::string &output) { try { this->connect(input, comp.getOutput(output)); } catch (LibMultiScaleException &e) { if (output == "") inputs[input] = comp; else throw e; } } /* --------------------------------------------------------------------- */ OutputContainer &Component::createOutput(const std::string &output) { outputs[output] = OutputContainer(*this, output); return outputs[output]; } /* --------------------------------------------------------------------- */ InputContainer &Component::createInput(const std::string &input) { inputs[input] = InputContainer(); return inputs[input]; } /* --------------------------------------------------------------------- */ void Component::removeInput(const std::string &input) { inputs.erase(input); } /* --------------------------------------------------------------------- */ void Component::removeOutput(const std::string &output) { outputs.erase(output); } /* --------------------------------------------------------------------- */ -void Component::compute() { +void Component::compute(bool allow_unconnected_failure) { auto secure_call = [&](auto &&f) { try { f(); + } catch (Component::UnconnectedInput &e) { + if (allow_unconnected_failure) { + DUMP("Unconnected input for component " << this->getID() << std::endl + << e.what(), + DBG_WARNING); + } else + throw e; + } catch (Component::NotInCommGroup &e) { + DUMP("This processor cannot not do anything as one input is generated " + "by another group of processors", + DBG_WARNING); } catch (LibMultiScaleException &e) { LM_FATAL_RE(e, "compute of " << this->getID() << " failed"); } catch (std::exception &e) { LM_FATAL("compute of " << this->getID() << " failed => " << e.what()); } catch (...) { LM_FATAL("compute of " << this->getID() << " failed => unknown exception was raised"); } }; bool need_recompute = this->evalDependencies(); bool amIinGroup = comm_group == nullptr || comm_group->amIinGroup(); if (not amIinGroup) { DUMP(this->getID() << ": Do not perform compute as I am not in the right group", DBG_INFO); throw NotInCommGroup{"For component '" + this->getID() + " I am not in the right group"}; calculated_once = true; return; } if (!need_recompute && calculated_once) return; // evalutate the inputs secure(this->evalInputs()); // get the context (comm_group, release) from inputs secure(this->acquireInputsContext()); // call the component code secure(this->compute_make_call()); // propagate the contexts from inputs to outputs secure(this->acquireInputsContext()); // sets the context to outputs secure(this->propagateContextToOutputs()); calculated_once = true; } #undef secure /* -------------------------------------------------------------------------- */ void Component::evalInputs() { for (auto &&[name, inp] : this->inputs) { if (not inp.has_value()) - LM_FATAL("Component(" << this->getID() << "): input(" << name - << ") was not created/computed/connected"); + throw Component::UnconnectedInput{"Component(" + this->getID() + + "): input(" + name + + ") was not created/computed/connected"}; + // LM_FATAL("Component(" << this->getID() << "): input(" << name + // << ") was not created/computed/connected"); inp.eval(); } } /* -------------------------------------------------------------------------- */ void Component::acquireInputsContext() { if (this->inputs.size() == 0) { this->changeRelease(); } for (auto &&[name, inp] : this->inputs) { if (not inp.has_value()) LM_FATAL(name << ": input was not created/computed/connected"); auto &connected_output = inp.get(); if (connected_output.has_argument()) this->acquireContext(connected_output.get()); } UInt release = this->evalRelease(); this->setRelease(release); } /* -------------------------------------------------------------------------- */ void Component::changeRelease() { LMObject::changeRelease(); for (auto &&[name, out] : this->outputs) { if (not out.has_argument()) continue; out->acquireContext(*this); } } /* -------------------------------------------------------------------------- */ bool Component::evalDependencies() { UInt release = this->evalRelease(); if (release != this->getRelease()) return true; return false; } /* -------------------------------------------------------------------------- */ void Component::propagateContextToOutputs() { for (auto &&[name, out] : this->outputs) { if (not out.has_value()) DUMP(name << ": output was not created/computed/connected", DBG_WARNING); try { out.get().acquireContext(*this); } catch (...) { DUMP(name << ": output could not receive context", DBG_WARNING); } } } /* -------------------------------------------------------------------------- */ UInt Component::evalRelease() { if (this->inputs.size() == 0) { return this->getRelease(); } std::vector release_inputs; for (auto &&[name, inp] : this->inputs) { UInt release; if (not inp.has_value()) return LMObject::random_release(); try { release = inp.evalRelease(); } catch (AutoDispatch::no_argument &e) { return LMObject::random_release(); } release_inputs.push_back(release); } UInt release = hash_vector(release_inputs); return release; } /* -------------------------------------------------------------------------- */ void Component::clear() { for (auto &&[name, out] : this->outputs) { this->applyArrayOutput([](auto &output) -> void { output.clear(); }, name); } } /* --------------------------------------------------------------------- */ UInt Component::size(const std::string &name) { try { auto &_output = this->getOutput(name); return _output.size(); } catch (AutoDispatch::no_argument &e) { // do nothing: cannot clean unallocated pointer } catch (AutoDispatch::bad_argument_cast &e) { // do nothing: it is not an array } catch (std::bad_cast &e) { // do nothing: it is not an array } return 0; } /* --------------------------------------------------------------------- */ __END_LIBMULTISCALE__ /* ---------------------------------------------------------------------- */ diff --git a/src/common/component_libmultiscale.hh b/src/common/component_libmultiscale.hh index 7ea6143..953a9cf 100644 --- a/src/common/component_libmultiscale.hh +++ b/src/common/component_libmultiscale.hh @@ -1,304 +1,304 @@ /** * @file component_libmultiscale.hh * * @author Guillaume Anciaux * * @date Mon Sep 08 23:40:22 2014 * * @brief This describe the root objects to be combined into valid LM * components * * @section LICENSE * * Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * LibMultiScale is free software: you can redistribute it and/or modify it * under the * terms of the GNU Lesser General Public License as published by the Free * Software Foundation, either version 3 of the License, or (at your option) any * later version. * * LibMultiScale is distributed in the hope that it will be useful, but * WITHOUT ANY * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR * A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more * details. * * You should have received a copy of the GNU Lesser General Public License * along with LibMultiScale. If not, see . * */ #ifndef __LIBMULTISCALE_COMPONENT_LIBMULTISCALE_HH__ #define __LIBMULTISCALE_COMPONENT_LIBMULTISCALE_HH__ /* -------------------------------------------------------------------------- */ #include "auto_arguments.hh" #include "container.hh" #include "lm_common.hh" #include "lm_object.hh" /* -------------------------------------------------------------------------- */ #include #include #include #include #include #include #include #include /* -------------------------------------------------------------------------- */ __BEGIN_LIBMULTISCALE__ /* -------------------------------------------------------------------------- */ class Component; /* -------------------------------------------------------------------------- */ /* -------------------------------------------------------------------------- */ // Output Container /* -------------------------------------------------------------------------- */ class OutputContainer : public AutoDispatch::ArgumentAny { public: OutputContainer() : component(nullptr){}; OutputContainer(Component &component, const std::string &name) : component(&component), name(name){}; OutputContainer(OutputContainer &&other) = default; OutputContainer(const OutputContainer &other) = default; inline OutputContainer &operator=(OutputContainer &&arg) = default; inline OutputContainer &operator=(OutputContainer &arg) = default; template decltype(auto) get() { return AutoDispatch::ArgumentAny::get(); } template decltype(auto) get() const { return AutoDispatch::ArgumentAny::get(); } using AutoDispatch::ArgumentAny::is_type; inline LMObject *operator->() { return &this->get(); } inline LMObject &operator*() { return this->get(); } inline const LMObject *operator->() const { return &this->get(); } inline const LMObject &operator*() const { return this->get(); } OutputContainer &operator=(Component &v); template inline OutputContainer &operator=(T &&v); template decltype(auto) alloc(T &&... construction_parameters); friend class InputContainer; UInt getRelease() const; UInt evalRelease() const; protected: Component *component; std::string name; }; /* --------------------------------------------------------------------- */ // InputContainer /* --------------------------------------------------------------------- */ class InputContainer { public: inline InputContainer() = default; inline InputContainer(InputContainer &&v) = default; inline InputContainer(const InputContainer &v) = default; inline InputContainer &operator=(InputContainer &&v) = default; InputContainer &operator=(OutputContainer &v); InputContainer &operator=(OutputContainer &&v); template inline InputContainer &operator=(T &&v); inline bool has_value() const { return output_reference.has_value(); } UInt getRelease() const; UInt evalRelease() const; OutputContainer &get() const; OutputContainer &eval() const; template inline decltype(auto) get(); template inline decltype(auto) get() const; protected: private: std::any output_reference; }; struct InputConnection { inline InputConnection(const std::string name) : name(name) {} template inline InputConnection &operator=(T &&t) { val = std::forward(t); return *this; } std::string name; OutputContainer val; }; inline InputConnection operator"" _input(const char *name, size_t) { return InputConnection(name); } /* --------------------------------------------------------------------- */ // Component /* --------------------------------------------------------------------- */ class Component : public virtual LMObject, public std::enable_shared_from_this { public: Component(); virtual ~Component(){}; //////////////////////////////////////////////////////////////// // meta-data of component //////////////////////////////////////////////////////////////// public: template void setCommGroup(CommGroup &group); void printself(std::ostream &os) const override; friend UInt _parse(Component &comp, std::stringstream &line, UInt); friend class OutputContainer; friend class InputContainer; //! increment the release void changeRelease() override; protected: struct UnconnectedInput : public std::runtime_error { using std::runtime_error::runtime_error; }; struct NotInCommGroup : public std::runtime_error { using std::runtime_error::runtime_error; }; //! evaluate all the inputs void evalInputs(); //! acquire the context fom all the inputs void acquireInputsContext(); //! returns true if the dependencies are recent => need recompute bool evalDependencies(); //! propagate the context to all the outputs void propagateContextToOutputs(); //! evaluate the current release (consult all chain of inputs) UInt evalRelease(); //////////////////////////////////////////////////////////////// // input management //////////////////////////////////////////////////////////////// public: template void connect(const std::string &input, T &&arg); void connect(const std::string &input, Component &comp, const std::string &output = ""); protected: inline const auto &getInputs() const { return inputs; } void clear_inputs() { inputs.clear(); } OutputContainer &getInput(const std::string &requested_input = ""); InputContainer &createInput(const std::string &input); private: InputContainer &getInternalInput(const std::string &requested_input = ""); //////////////////////////////////////////////////////////////// // output management //////////////////////////////////////////////////////////////// public: template inline decltype(auto) evalArrayOutput(const std::string &output_name = ""); template inline decltype(auto) evalOutput(const std::string &output_name = ""); OutputContainer &evalOutput(const std::string &output_name = ""); const std::map &evalOutputs(); protected: const std::map &getOutputs(); template ContainerArray &getOutputAsArray(const std::string &output_name = ""); OutputContainer &createOutput(const std::string &output); template inline decltype(auto) getOutput(const std::string &output_name = ""); OutputContainer &getOutput(const std::string &output_name = ""); template void createArrayOutputs(Names &&output_names); template void createArrayOutput(const std::string &name); template auto &allocOutput(const std::string &output_name, T &&... construction_parameters); void removeInput(const std::string &input); void removeOutput(const std::string &output); //////////////////////////////////////////////////////////////// // compute management //////////////////////////////////////////////////////////////// public: virtual void compute_make_call() = 0; - virtual void compute(); + virtual void compute(bool allow_unconnected_failure = false); template std::enable_if_t, InputConnection>> compute(T &&arg, Ts &&... args); template std::enable_if_t, InputConnection>> compute(T &&arg); private: template decltype(auto) getContainer(const std::string &name, ContMap &cont); //////////////////////////////////////////////////////////////// // array output management //////////////////////////////////////////////////////////////// public: void clear(); UInt size(const std::string &name = ""); private: template inline void applyArrayOutput(F &&f, const std::string &name = ""); protected: std::map inputs; std::map outputs; std::string default_output = ""; bool calculated_once; }; /* --------------------------------------------------------------------- */ template class SingleOutputComponent : public Component { public: using Component::Component; SingleOutputComponent(const LMID &name) : LMObject(name), Component(name) {} operator Output &() { return this->getOutputs().begin()->second.get().template cast(); } auto begin() { return (*this)->begin(); } auto end() { return (*this)->end(); } Output *operator->() { return &static_cast(*this); } Output &operator*() { return *this; } }; template struct casted_component { template using cast = SingleOutputComponent; }; /* --------------------------------------------------------------------- */ __END_LIBMULTISCALE__ /* --------------------------------------------------------------------- */ #include "component_libmultiscale_inline_impl.hh" #include "possible_types.hh" /* --------------------------------------------------------------------- */ #endif /* __LIBMULTISCALE_COMPONENT_LIBMULTISCALE_HH__ */ diff --git a/src/coupling/xiao.cc b/src/coupling/xiao.cc index 201b0ad..2c8f88c 100644 --- a/src/coupling/xiao.cc +++ b/src/coupling/xiao.cc @@ -1,345 +1,349 @@ /** * @file xiao.cc * * @author Guillaume Anciaux * * @date Wed Jan 15 17:00:43 2014 * * @brief Bridging Domain/Arlequin method * * @section LICENSE * * Copyright INRIA and CEA * * The LibMultiScale is a C++ parallel framework for the multiscale * coupling methods dedicated to material simulations. This framework * provides an API which makes it possible to program coupled simulations * and integration of already existing codes. * * This Project was initiated in a collaboration between INRIA Futurs Bordeaux * within ScAlApplix team and CEA/DPTA Ile de France. * The project is now continued at the Ecole Polytechnique Fédérale de Lausanne * within the LSMS/ENAC laboratory. * * This software is governed by the CeCILL-C license under French law and * abiding by the rules of distribution of free software. You can use, * modify and/ or redistribute the software under the terms of the CeCILL-C * license as circulated by CEA, CNRS and INRIA at the following URL * "http://www.cecill.info". * * As a counterpart to the access to the source code and rights to copy, * modify and redistribute granted by the license, users are provided only * with a limited warranty and the software's author, the holder of the * economic rights, and the successive licensors have only limited * liability. * * In this respect, the user's attention is drawn to the risks associated * with loading, using, modifying and/or developing or reproducing the * software by the user in light of its specific status of free software, * that may mean that it is complicated to manipulate, and that also * therefore means that it is reserved for developers and experienced * professionals having in-depth computer knowledge. Users are therefore * encouraged to load and test the software's suitability as regards their * requirements in conditions enabling the security of their systems and/or * data to be ensured and, more generally, to use and operate it in the * same conditions as regards security. * * The fact that you are presently reading this means that you have had * knowledge of the CeCILL-C license and that you accept its terms. * */ #define TIMER /* -------------------------------------------------------------------------- */ #include "xiao.hh" #include "lib_bridging.hh" #include "lm_common.hh" #include "stimulation_zero.hh" /* -------------------------------------------------------------------------- */ __BEGIN_LIBMULTISCALE__ /* -------------------------------------------------------------------------- */ Xiao::Xiao(const std::string &name) : LMObject(name), ArlequinTemplate(name) { this->createOutput("unmatched-point-boundary"); this->createOutput("unmatched-point-bridging"); this->createOutput("unmatched-mesh-boundary"); this->createOutput("unmatched-mesh-bridging"); this->createOutput("matched-point-boundary"); this->createOutput("matched-point-bridging"); this->createOutput("matched-mesh-boundary"); this->createOutput("matched-mesh-bridging"); this->createOutput("rhs"); this->createOutput("A"); } /* -------------------------------------------------------------------------- */ Xiao::~Xiao() {} /* -------------------------------------------------------------------------- */ void Xiao::coupling(CouplingStage stage) { if (!this->is_in_continuum() && !this->is_in_atomic()) return; if (stage == COUPLING_STEP2) { DUMP("Zero boundary condition", DBG_INFO); STARTTIMER("ZERO_BC"); boundary_zone.setPointField(_force, 0., spatial_dimension); STOPTIMER("ZERO_BC"); } if (stage != COUPLING_STEP3) return; STARTTIMER("BuildRHS"); DUMP("Build RHS", DBG_INFO); this->buildRHS(_velocity); STOPTIMER("BuildRHS"); DUMP("Synch with migration bridging", DBG_INFO); STARTTIMER("updatesForMigrations"); bridging_zone.updateForMigration(); // this->size_constraint = bridging_zone.getNumberLocalMatchedPoints(); DUMP("Synch with migration boundary", DBG_INFO); boundary_zone.updateForMigration(); STOPTIMER("updatesForMigrations"); DUMP("Correcting surface effect", DBG_INFO); STARTTIMER("projectAtomicVelocitiesOnMe.hh"); // Parent::MDboundary_zone.projectAtomicDOFsOnMesh(_displacement); boundary_zone.projectAtomicFieldOnMesh(_velocity); STOPTIMER("projectAtomicVelocitiesOnMe.hh"); STARTTIMER("BuildRHS"); DUMP("Build RHS", DBG_INFO); this->buildRHS(_velocity); STOPTIMER("BuildRHS"); DUMP("solving constraint", DBG_INFO); STARTTIMER("Solving constraint"); this->solveConstraint(); STOPTIMER("Solving constraint"); DUMP("Applying constraint", DBG_INFO); STARTTIMER("Correcting"); this->applyCorrection(_velocity); STOPTIMER("Correcting"); DUMP(this->getID() << " : coupling stage all done", DBG_INFO); } /* -------------------------------------------------------------------------- */ void Xiao::solveConstraint() { for (auto &&[r, a] : zip(rhs.rowwise(), A)) { r /= a; } } /* -------------------------------------------------------------------------- */ void Xiao::applyCorrection(FieldType field) { if (this->is_in_continuum()) { auto &correction = this->bridging_zone.buffer_for_nodes; correction.resize(rhs.rows(), rhs.cols()); correction = -(this->bridging_zone.smatrix * this->rhs.matrix()).array(); for (auto &&[corr, lbda] : zip(correction.rowwise(), this->lambdas_mesh)) { corr /= lbda; } this->bridging_zone.addMeshField(field, correction); } if (this->is_in_atomic()) { auto &correction = this->bridging_zone.buffer_for_points; correction.resize(rhs.rows(), rhs.cols()); for (auto &&[cor, r, lbda] : zip(correction.rowwise(), rhs.rowwise(), lambdas_point)) { LM_ASSERT(lbda, "weight associated with atom is zero : abort"); cor = r / lbda; } this->bridging_zone.addPointField(field, correction); } } /* -------------------------------------------------------------------------- */ void Xiao::buildRHS(FieldType field_type) { DUMP("Clean RHS", DBG_INFO); if (this->rhs.rows() == 0) this->rhs.resize(0, spatial_dimension); this->rhs.setZero(); if (this->is_in_continuum()) { this->bridging_zone.mesh2Point(field_type, this->rhs); } if (this->is_in_atomic()) { ComputeExtract field(this->getID()); field.setParam("FIELD", field_type); field.compute(this->bridging_zone.pointList); - this->rhs -= field.evalArrayOutput().array(); + if (not this->is_in_continuum()) { + this->rhs = field.evalArrayOutput().array(); + } else { + this->rhs -= field.evalArrayOutput().array(); + } } bridging_zone.synchronizeVectorBySum("rhs", this->rhs); } /* -------------------------------------------------------------------------- */ template void Xiao::init(DomainA &domA, DomainC &domC) { ArlequinTemplate::init(domA, domC); // allocate the bridging zone (parallel version) bridging_zone.setParam("GEOMETRY", this->bridging_geom); bridging_zone.setParam("CHECK_COHERENCY", this->check_coherency); bridging_zone.setParam("CENTROID", this->centroid_flag); // set the pbc pairs bridging_zone.setPBCPairs(domC.getPBCpairs()); // initialize the bridging_zone object bridging_zone.init(domA, domC); // allocate the vectors necessary to continue if (this->is_in_atomic() or bridging_zone.getNumberLocalMatchedPoints()) { this->bridging_zone.attachVector(A); // this->allocate(bridging_zone.getNumberLocalMatchedPoints()); } // build weights this->computeWeights(bridging_zone.pointList, bridging_zone.meshList); /* build constraint matrix */ this->buildConstraintMatrix(); /* now treat the boundary zone */ boundary_zone.setParam("GEOMETRY", this->boundary_geom); boundary_zone.setParam("CHECK_COHERENCY", this->check_coherency); boundary_zone.setParam("CENTROID", this->centroid_flag); // initialize the bridging_zone object boundary_zone.init(domA, domC); // this->boundary_zone.projectAtomicDOFsOnMesh(_displacement); boundary_zone.projectAtomicFieldOnMesh(_velocity); boundary_zone.setPointField(_force, 0., spatial_dimension); this->getOutput("unmatched-point-boundary") = boundary_zone.unmatchedPointList; this->getOutput("unmatched-point-bridging") = bridging_zone.unmatchedPointList; this->getOutput("unmatched-mesh-boundary") = boundary_zone.unmatchedMeshList; this->getOutput("unmatched-mesh-bridging") = bridging_zone.unmatchedMeshList; this->getOutput("matched-point-boundary") = boundary_zone.pointList; this->getOutput("matched-point-bridging") = bridging_zone.pointList; this->getOutput("matched-mesh-boundary") = boundary_zone.meshList; this->getOutput("matched-mesh-bridging") = bridging_zone.meshList; this->getOutput("rhs") = rhs; this->getOutput("A") = A; } /* -------------------------------------------------------------------------- */ void Xiao::buildConstraintMatrix() { if (this->is_in_continuum()) { auto &shp = this->bridging_zone.smatrix; Array lumped_shape(shp.rows(), 1); for (UInt i = 0; i < shp.rows(); ++i) { lumped_shape(i) = shp.row(i).sum(); } A = (shp.transpose() * (1. / lambdas_mesh * lumped_shape).matrix()); } if (this->is_in_atomic()) { if (not this->is_in_continuum()) A = Array::Zero(lambdas_point.rows(), lambdas_point.cols()); A += 1. / lambdas_point; } /* synch constraint matrix */ bridging_zone.synchronizeVectorBySum("A", this->A); } /* -------------------------------------------------------------------------- */ /* LMDESC XIAO This class implements the BridgingDomain/Arlequin method. Two zones are declared, one for the coupling (light blue atoms) one for providing a stiff boundary condition to the atoms (purple atoms). It defines a coupling area as depicted below: .. image:: images/bridging-zone-BM.svg In the bridging part a linear weight function is built the weight the energies in the Arlequin method flavor. The detailed algorithm is: .. math:: \left\{ \begin{array}{l} \displaystyle\hat{M}_I \ddot{\mathbf{u}}_I = - \mathbf{f}^R_I + \sum_{k=1}^L \mathbf{\lambda}_k \frac{\partial \mathbf{g}_k}{\partial \mathbf{u}_I} \\ \displaystyle \hat{m}_i \ddot{\mathbf{d}}_i = \mathbf{f}^R_i + \sum_{k=1}^L \mathbf{\lambda}_k \frac{\partial \mathbf{g}_k}{\partial \mathbf{d}_i} \end{array} \right. where :math:`\alpha(X_I)` and :math:`\alpha(X_i)` are the weight associated with nodes and atoms respectively, where :math:`\hat{M}_I=\alpha(X_I)M_I` and :math:`\hat{m}_i=(1-\alpha(X_i))m_i` and where the constraint multipliers are computed by solving :math:`H \Lambda = \mathbf{g}^\star` with: .. math:: H_{ik} = \Delta t \left( \sum_J \varphi_J(\mathbf{X}_I) \hat{M_J}^{-1} \frac{\partial \mathbf{g}_k}{\partial \mathbf{u}_J} - \hat{m}_I^{-1} \frac{\partial \mathbf{g}_k}{\partial \mathbf{d}_i} \right), \mathbf{g}_i^\star = \sum_J \varphi_J(\mathbf{X}_I) \dot{\mathbf{u}}_J^\star - \dot{\mathbf{d}}_i^\star, where :math:`u^\star, d^\star` are the displacement obtained after one time step by the integration scheme when the constraints are not applied. For debugging puposes several computes are registered to the central system: - weight-fe-\${couplerID} : it is a compute containing the weights associated with nodes inside of the bridging zone (the :math:`\alpha(X_I)`). - weight-md-\${couplerID} : it is a compute containing the weights associated with atoms inside of the bridging zone (the :math:`\alpha(X_i)`). - lambdas-fe-\${couplerID} : it is a compute containing the :math:`\hat{M}_I`. - lambdas-md-\${couplerID} : it is a compute containing the :math:`\hat{m}_i`. */ /* LMEXAMPLE COUPLING_CODE bdmID md fe XIAO GEOMETRY 3 BOUNDARY 4 GRID_DIVISIONX 10 GRID_DIVISIONY 10 QUALITY 1e-3 */ /* LMHERITANCE arlequin_template */ void Xiao::declareParams() { ArlequinTemplate::declareParams(); /* LMKEYWORD CENTROID Set the selection of the bridging zones based on centroid. */ this->parseTag("CENTROID", centroid_flag, false); } /* -------------------------------------------------------------------------- */ DECLARE_COUPLER_INIT_MAKE_CALL(Xiao, domA, domC) /* -------------------------------------------------------------------------- */ __END_LIBMULTISCALE__ diff --git a/src/dumper/dumper_interface.hh b/src/dumper/dumper_interface.hh index 4813abb..14884b9 100644 --- a/src/dumper/dumper_interface.hh +++ b/src/dumper/dumper_interface.hh @@ -1,133 +1,123 @@ /** * @file dumper_interface.hh * * @author Guillaume Anciaux * * @date Wed Apr 03 23:50:25 2013 * * @brief The interface shared by all dumpers * * @section LICENSE * * Copyright INRIA and CEA * * The LibMultiScale is a C++ parallel framework for the multiscale * coupling methods dedicated to material simulations. This framework * provides an API which makes it possible to program coupled simulations * and integration of already existing codes. * * This Project was initiated in a collaboration between INRIA Futurs Bordeaux * within ScAlApplix team and CEA/DPTA Ile de France. * The project is now continued at the Ecole Polytechnique Fédérale de Lausanne * within the LSMS/ENAC laboratory. * * This software is governed by the CeCILL-C license under French law and * abiding by the rules of distribution of free software. You can use, * modify and/ or redistribute the software under the terms of the CeCILL-C * license as circulated by CEA, CNRS and INRIA at the following URL * "http://www.cecill.info". * * As a counterpart to the access to the source code and rights to copy, * modify and redistribute granted by the license, users are provided only * with a limited warranty and the software's author, the holder of the * economic rights, and the successive licensors have only limited * liability. * * In this respect, the user's attention is drawn to the risks associated * with loading, using, modifying and/or developing or reproducing the * software by the user in light of its specific status of free software, * that may mean that it is complicated to manipulate, and that also * therefore means that it is reserved for developers and experienced * professionals having in-depth computer knowledge. Users are therefore * encouraged to load and test the software's suitability as regards their * requirements in conditions enabling the security of their systems and/or * data to be ensured and, more generally, to use and operate it in the * same conditions as regards security. * * The fact that you are presently reading this means that you have had * knowledge of the CeCILL-C license and that you accept its terms. * */ #ifndef __LIBMULTISCALE_DUMPER_INTERFACE_HH__ #define __LIBMULTISCALE_DUMPER_INTERFACE_HH__ /* -------------------------------------------------------------------------- */ #include "action_interface.hh" extern std::string lm_release_info; /* -------------------------------------------------------------------------- */ __BEGIN_LIBMULTISCALE__ /** Class DumperInterface. * Interface class to all dumpers */ class DumperInterface : public ActionInterface { /* ------------------------------------------------------------------------ */ /* Constructors/Destructors */ /* ------------------------------------------------------------------------ */ public: DumperInterface(); virtual ~DumperInterface(){}; /* ------------------------------------------------------------------------ */ /* Methods */ /* ------------------------------------------------------------------------ */ //! set parameters virtual void declareParams(); //! construct the name of the files to be generated void buildNameFile(); //! init stuff virtual void init(); const std::string &getBaseName(); //! call for immediate dump void dump(); /* ------------------------------------------------------------------------ */ /* Class Members */ /* ------------------------------------------------------------------------ */ protected: //! flag to set if no action required from dumper bool nothing_to_do_flag; //! data stored by step std::map> by_step_data; //! unit system in which to output the dump data, not used by all dumpers UnitSystem dump_unit_system; //! domain id to apply stimulator on UInt domain_id; //! prefix for directory std::string prefix; //! base name of the filename std::string base_name; }; /* -------------------------------------------------------------------------- */ #define DECLARE_DUMPER(class_name, ...) \ class_name(const std::string &name); \ virtual ~class_name(); \ DECORATE_FUNCTION_DISPATCH(dump, __VA_ARGS__) \ EIGEN_MAKE_ALIGNED_OPERATOR_NEW \ \ void declareParams() override; \ void compute_make_call() override #define DECLARE_DUMPER_MAKE_CALL(class_name) \ void class_name::compute_make_call() { \ - try { \ - this->dump(this->getInput()); \ - } catch (Component::UnconnectedInput & e) { \ - LM_FATAL("Unconnected input for component " << this->getID() \ - << std::endl \ - << e.what()); \ - } catch (Component::NotInCommGroup & e) { \ - DUMP("This processor cannot not do anything as one input is generated " \ - "by another group of processors", \ - DBG_WARNING); \ - } \ + this->dump(this->getInput()); \ } __END_LIBMULTISCALE__ #endif /* __LIBMULTISCALE_DUMPER_INTERFACE_HH__ */ diff --git a/src/filter/filter_interface.hh b/src/filter/filter_interface.hh index 5770667..0eb9d24 100644 --- a/src/filter/filter_interface.hh +++ b/src/filter/filter_interface.hh @@ -1,92 +1,86 @@ /** * @file filter_interface.hh * * * @date Mon Oct 28 19:23:14 2013 * * @brief This is the interface to all filters * * @section LICENSE * * Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * LibMultiScale is free software: you can redistribute it and/or modify it * under the * terms of the GNU Lesser General Public License as published by the Free * Software Foundation, either version 3 of the License, or (at your option) any * later version. * * LibMultiScale is distributed in the hope that it will be useful, but * WITHOUT ANY * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR * A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more * details. * * You should have received a copy of the GNU Lesser General Public License * along with LibMultiScale. If not, see . * */ /* ---------------------------------------------------------------- */ #ifndef __LIBMULTISCALE_FILTER_INTERFACE_HH__ #define __LIBMULTISCALE_FILTER_INTERFACE_HH__ /* ---------------------------------------------------------------- */ #include "action_interface.hh" #include "component_libmultiscale.hh" /* ---------------------------------------------------------------- */ __BEGIN_LIBMULTISCALE__ /* ---------------------------------------------------------------- */ class ActionInterface; template class FilterCompatibility; /* ---------------------------------------------------------------- */ class FilterInterface : public ActionInterface { public: FilterInterface(){}; virtual ~FilterInterface(){}; virtual void declareParams() override; protected: DOFType dt; }; /* ---------------------------------------------------------------- */ /// standard output stream operator inline std::ostream &operator<<(std::ostream &stream, const FilterInterface &_this) { _this.printself(stream); return stream; } __END_LIBMULTISCALE__ /* ---------------------------------------------------------------- */ #define DECLARE_FILTER(class_name, ...) \ class_name(const std::string &name); \ virtual ~class_name(); \ \ private: \ DECORATE_FUNCTION_DISPATCH(build, __VA_ARGS__) \ EIGEN_MAKE_ALIGNED_OPERATOR_NEW \ \ public: \ void declareParams() override; \ void compute_make_call() override #define DECLARE_COMPUTE_MAKE_CALL(class_name) \ void class_name::compute_make_call() { \ - try { \ - DUMP(this->getID() << ": Compute", DBG_INFO); \ - this->build(this->getInput()); \ - } catch (Component::UnconnectedInput & e) { \ - LM_FATAL("Unconnected input for component " << this->getID() \ - << std::endl \ - << e.what()); \ - } \ + DUMP(this->getID() << ": Compute", DBG_INFO); \ + this->build(this->getInput()); \ } #endif /* __LIBMULTISCALE_FILTER_INTERFACE_HH__ */ diff --git a/src/stimulation/stimulation_interface.hh b/src/stimulation/stimulation_interface.hh index 41d53d9..fe43891 100644 --- a/src/stimulation/stimulation_interface.hh +++ b/src/stimulation/stimulation_interface.hh @@ -1,112 +1,102 @@ /** * @file stimulation_interface.hh * * @author Guillaume Anciaux * * @date Wed Apr 03 23:50:25 2013 * * @brief This is the interface of all stimulators * * @section LICENSE * * Copyright INRIA and CEA * * The LibMultiScale is a C++ parallel framework for the multiscale * coupling methods dedicated to material simulations. This framework * provides an API which makes it possible to program coupled simulations * and integration of already existing codes. * * This Project was initiated in a collaboration between INRIA Futurs Bordeaux * within ScAlApplix team and CEA/DPTA Ile de France. * The project is now continued at the Ecole Polytechnique Fédérale de Lausanne * within the LSMS/ENAC laboratory. * * This software is governed by the CeCILL-C license under French law and * abiding by the rules of distribution of free software. You can use, * modify and/ or redistribute the software under the terms of the CeCILL-C * license as circulated by CEA, CNRS and INRIA at the following URL * "http://www.cecill.info". * * As a counterpart to the access to the source code and rights to copy, * modify and redistribute granted by the license, users are provided only * with a limited warranty and the software's author, the holder of the * economic rights, and the successive licensors have only limited * liability. * * In this respect, the user's attention is drawn to the risks associated * with loading, using, modifying and/or developing or reproducing the * software by the user in light of its specific status of free software, * that may mean that it is complicated to manipulate, and that also * therefore means that it is reserved for developers and experienced * professionals having in-depth computer knowledge. Users are therefore * encouraged to load and test the software's suitability as regards their * requirements in conditions enabling the security of their systems and/or * data to be ensured and, more generally, to use and operate it in the * same conditions as regards security. * * The fact that you are presently reading this means that you have had * knowledge of the CeCILL-C license and that you accept its terms. * */ /* -------------------------------------------------------------------------- */ #ifndef __LIBMULTISCALE_STIMULATION_INTERFACE_HH__ #define __LIBMULTISCALE_STIMULATION_INTERFACE_HH__ /* -------------------------------------------------------------------------- */ #include "action_interface.hh" /* -------------------------------------------------------------------------- */ __BEGIN_LIBMULTISCALE__ class StimulationInterface : public ActionInterface { public: StimulationInterface(); virtual ~StimulationInterface(); void init() override{}; void stimulate(); protected: }; /* -------------------------------------------------------------------------- */ #define DECLARE_STIMULATION(class_name, ...) \ class_name(const std::string &name); \ virtual ~class_name(); \ DECORATE_FUNCTION_DISPATCH(stimulate, __VA_ARGS__) \ EIGEN_MAKE_ALIGNED_OPERATOR_NEW \ \ void declareParams() override; \ void compute_make_call() override decltype(auto) inline make_input_names() { return std::make_tuple(""); } template decltype(auto) make_input_names(Args... args) { return std::make_tuple(args...); } /* -------------------------------------------------------------------------- */ #define DECLARE_STIMULATION_MAKE_CALL(class_name, ...) \ void class_name::compute_make_call() { \ auto input_names = make_input_names(stringify_macro(__VA_ARGS__)); \ - try { \ - std::apply( \ - [&](auto &&... input) { \ - this->stimulate(this->getInput(input)...); \ - }, \ - input_names); \ - } catch (Component::UnconnectedInput & e) { \ - LM_FATAL("Unconnected input for component " << this->getID() \ - << std::endl \ - << e.what()); \ - } catch (Component::NotInCommGroup & e) { \ - DUMP("This processor cannot not do anything as one input is generated " \ - "by another group of processors", \ - DBG_WARNING); \ - } \ + std::apply( \ + [&](auto &&... input) { \ + this->stimulate(this->getInput(input)...); \ + }, \ + input_names); \ } __END_LIBMULTISCALE__ #endif /* __LIBMULTISCALE_STIMULATION_INTERFACE_HH__ */