source: src/Actions/MoleculeAction/ForceAnnealingAction.cpp@ 322d58

AutomationFragmentation_failures Candidate_v1.6.1 ChemicalSpaceEvaluator Enhanced_StructuralOptimization Enhanced_StructuralOptimization_continued Exclude_Hydrogens_annealWithBondGraph ForceAnnealing_with_BondGraph ForceAnnealing_with_BondGraph_contraction-expansion Gui_displays_atomic_force_velocity JobMarket_RobustOnKillsSegFaults JobMarket_StableWorkerPool PythonUI_with_named_parameters StoppableMakroAction TremoloParser_IncreasedPrecision
Last change on this file since 322d58 was 322d58, checked in by Frederik Heber <frederik.heber@…>, 7 years ago

ForceAnnealing now uses step width according to Barzilai-Borwein method.

  • this is not as good as what MPQC does but as long as Actions don't have an internal state, we cannot do anything better, e.g. a real line-search.
  • removed paramater delta from optimize-structure.
  • DOCU: Updated entry in user guide and noted Barzilai-Borwein method.
  • TESTS: removed deltat from force annealing regression tests.
  • TESTS: switched off StructureOptimization integration tests as poolworkers crash at the moment, making the tests hang indefinitely.
  • Property mode set to 100644
File size: 5.1 KB
Line 
1/*
2 * Project: MoleCuilder
3 * Description: creates and alters molecular systems
4 * Copyright (C) 2014 Frederik Heber. All rights reserved.
5 *
6 *
7 * This file is part of MoleCuilder.
8 *
9 * MoleCuilder is free software: you can redistribute it and/or modify
10 * it under the terms of the GNU General Public License as published by
11 * the Free Software Foundation, either version 2 of the License, or
12 * (at your option) any later version.
13 *
14 * MoleCuilder is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 * GNU General Public License for more details.
18 *
19 * You should have received a copy of the GNU General Public License
20 * along with MoleCuilder. If not, see <http://www.gnu.org/licenses/>.
21 */
22
23/*
24 * ForceAnnealingAction.cpp
25 *
26 * Created on: Aug 02, 2014
27 * Author: heber
28 */
29
30// include config.h
31#ifdef HAVE_CONFIG_H
32#include <config.h>
33#endif
34
35//#include "CodePatterns/MemDebug.hpp"
36
37#include "Actions/UndoRedoHelpers.hpp"
38#include "Atom/atom.hpp"
39#include "Atom/AtomicInfo.hpp"
40#include "Atom/AtomSet.hpp"
41#include "CodePatterns/Log.hpp"
42#include "CodePatterns/Verbose.hpp"
43#include "Dynamics/ForceAnnealing.hpp"
44#include "molecule.hpp"
45#include "World.hpp"
46#include "WorldTime.hpp"
47
48#include <vector>
49#include <iostream>
50#include <fstream>
51#include <string>
52
53#include "Actions/MoleculeAction/ForceAnnealingAction.hpp"
54
55using namespace MoleCuilder;
56
57enum VectorIndexType {
58 PositionIndex=0,
59 VelocityIndex=1,
60 ForceIndex=2
61};
62
63// and construct the stuff
64#include "ForceAnnealingAction.def"
65#include "Action_impl_pre.hpp"
66/** =========== define the function ====================== */
67ActionState::ptr MoleculeForceAnnealingAction::performCall() {
68 AtomSetMixin<std::vector<atom *> > set(World::getInstance().getSelectedAtoms());
69 if (set.empty()) {
70 STATUS("No atoms selected.");
71 return Action::failure;
72 }
73 // we always operate relative to current time step, except on single debug output
74 size_t CurrentStep = WorldTime::getInstance().getTime();
75 if (params.DoOutput.get()) {
76 // copy current time step to new one and and proceed on this one
77 for (World::AtomSelectionIterator iter = World::getInstance().beginAtomSelection();
78 iter != World::getInstance().endAtomSelection();
79 ++iter) {
80 atom * const Walker = iter->second;
81 Walker->setPositionAtStep(CurrentStep+1,
82 Walker->getPositionAtStep(CurrentStep));
83 Walker->setAtomicVelocityAtStep(CurrentStep+1,
84 Walker->getAtomicVelocityAtStep(CurrentStep));
85 Walker->setAtomicForceAtStep(CurrentStep+1,
86 Walker->getAtomicForceAtStep(CurrentStep));
87 }
88 // increment to next time step: re-creates bond graph
89 ++CurrentStep;
90 World::getInstance().setTime(CurrentStep);
91 }
92 ForceAnnealing<std::vector<atom *> > optimizer(
93 set,
94 true,
95 params.steps.get());
96 // parse forces into next step
97 if (!params.forcesfile.get().string().empty()) {
98 LOG(1, "Parsing forces file.");
99 if (!optimizer.parseForcesFile(params.forcesfile.get().string().c_str(), CurrentStep))
100 LOG(2, "File " << params.forcesfile.get() << " not found.");
101 else
102 LOG(2, "File " << params.forcesfile.get() << " found and parsed.");
103 }
104
105 // create undo state for all selected atoms (undo info)
106 std::vector<AtomicInfo> UndoInfo;
107 UndoInfo.reserve(set.size());
108 {
109 for (World::AtomSelectionConstIterator iter = World::getInstance().beginAtomSelection();
110 iter != World::getInstance().endAtomSelection();
111 ++iter)
112 UndoInfo.push_back(AtomicInfo(*(iter->second)));
113 }
114
115 // perform optimization step
116 LOG(1, "Structural optimization.");
117 optimizer(CurrentStep, 1);
118 STATUS("Successfully optimized structure by one step.");
119
120// // increment to next time step
121// World::getInstance().setTime(CurrentStep+1);
122
123 std::vector<AtomicInfo> RedoInfo;
124 RedoInfo.reserve(set.size());
125 {
126 for (World::AtomSelectionConstIterator iter = World::getInstance().beginAtomSelection();
127 iter != World::getInstance().endAtomSelection();
128 ++iter)
129 RedoInfo.push_back(AtomicInfo(*(iter->second)));
130 }
131 MoleculeForceAnnealingState *UndoState =
132 new MoleculeForceAnnealingState(UndoInfo, RedoInfo, params);
133
134 return ActionState::ptr(UndoState);
135}
136
137ActionState::ptr MoleculeForceAnnealingAction::performUndo(ActionState::ptr _state) {
138 MoleculeForceAnnealingState *state =
139 assert_cast<MoleculeForceAnnealingState*>(_state.get());
140
141 // set stored old state
142 SetAtomsFromAtomicInfo(state->UndoInfo);
143
144 return ActionState::ptr(_state);
145}
146
147ActionState::ptr MoleculeForceAnnealingAction::performRedo(ActionState::ptr _state){
148 MoleculeForceAnnealingState *state =
149 assert_cast<MoleculeForceAnnealingState*>(_state.get());
150
151 // set stored new state
152 SetAtomsFromAtomicInfo(state->RedoInfo);
153
154 return ActionState::ptr(_state);
155}
156
157bool MoleculeForceAnnealingAction::canUndo() {
158 return true;
159}
160
161bool MoleculeForceAnnealingAction::shouldUndo() {
162 return true;
163}
164/** =========== end of function ====================== */
Note: See TracBrowser for help on using the repository browser.