[d2a0f6d] | 1 | /*
|
---|
| 2 | * Project: MoleCuilder
|
---|
| 3 | * Description: creates and alters molecular systems
|
---|
| 4 | * Copyright (C) 2012 University of Bonn. All rights reserved.
|
---|
[5aaa43] | 5 | * Copyright (C) 2013 Frederik Heber. All rights reserved.
|
---|
[d2a0f6d] | 6 | *
|
---|
| 7 | *
|
---|
| 8 | * This file is part of MoleCuilder.
|
---|
| 9 | *
|
---|
| 10 | * MoleCuilder is free software: you can redistribute it and/or modify
|
---|
| 11 | * it under the terms of the GNU General Public License as published by
|
---|
| 12 | * the Free Software Foundation, either version 2 of the License, or
|
---|
| 13 | * (at your option) any later version.
|
---|
| 14 | *
|
---|
| 15 | * MoleCuilder is distributed in the hope that it will be useful,
|
---|
| 16 | * but WITHOUT ANY WARRANTY; without even the implied warranty of
|
---|
| 17 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
---|
| 18 | * GNU General Public License for more details.
|
---|
| 19 | *
|
---|
| 20 | * You should have received a copy of the GNU General Public License
|
---|
| 21 | * along with MoleCuilder. If not, see <http://www.gnu.org/licenses/>.
|
---|
| 22 | */
|
---|
| 23 |
|
---|
| 24 | /*
|
---|
| 25 | * VMGJob.cpp
|
---|
| 26 | *
|
---|
| 27 | * Created on: Jul 13, 2012
|
---|
| 28 | * Author: heber
|
---|
| 29 | */
|
---|
| 30 |
|
---|
| 31 |
|
---|
| 32 | // include config.h
|
---|
| 33 | #ifdef HAVE_CONFIG_H
|
---|
| 34 | #include <config.h>
|
---|
| 35 | #endif
|
---|
| 36 |
|
---|
[8a8c8c] | 37 | #ifdef HAVE_MPI
|
---|
| 38 | #include "mpi.h"
|
---|
| 39 | #endif
|
---|
| 40 |
|
---|
[d2a0f6d] | 41 | // include headers that implement a archive in simple text format
|
---|
| 42 | // otherwise BOOST_CLASS_EXPORT_IMPLEMENT has no effect
|
---|
| 43 | #include <boost/archive/text_oarchive.hpp>
|
---|
| 44 | #include <boost/archive/text_iarchive.hpp>
|
---|
| 45 |
|
---|
[e9cfc4] | 46 | #include "mg.hpp"
|
---|
| 47 | #include "base/object.hpp"
|
---|
| 48 | #include "base/defs.hpp"
|
---|
[8a8c8c] | 49 |
|
---|
| 50 | #ifdef HAVE_MPI
|
---|
| 51 | #include "comm/comm_mpi.hpp"
|
---|
| 52 | #include "comm/domain_decomposition_mpi.hpp"
|
---|
| 53 | #else
|
---|
[e9cfc4] | 54 | #include "comm/comm_serial.hpp"
|
---|
[8a8c8c] | 55 | #endif
|
---|
[ee9018] | 56 | // periodic boundary conditions
|
---|
[e9cfc4] | 57 | #include "cycles/cycle_cs_periodic.hpp"
|
---|
[e8e472b] | 58 | #include "discretization/discretization_poisson_fd.hpp"
|
---|
[ee9018] | 59 | #include "level/level_operator_cs.hpp"
|
---|
| 60 | #include "smoother/gsrb_poisson_4.hpp"
|
---|
| 61 | #include "solver/solver_singular.hpp"
|
---|
| 62 | // open boundary conditions
|
---|
| 63 | #include "cycles/cycle_fas_dirichlet.hpp"
|
---|
| 64 | #include "discretization/discretization_poisson_fv.hpp"
|
---|
[e9cfc4] | 65 | #include "grid/multigrid.hpp"
|
---|
| 66 | //#include "grid/tempgrid.hpp"
|
---|
[ee9018] | 67 | #include "level/level_operator_fas.hpp"
|
---|
[e9cfc4] | 68 | #include "level/stencils.hpp"
|
---|
[ee9018] | 69 | #include "smoother/gsrb_poisson_2.hpp"
|
---|
[e9cfc4] | 70 | #include "solver/givens.hpp"
|
---|
[ee9018] | 71 | #include "solver/solver_regular.hpp"
|
---|
[a82602] | 72 | #include "units/particle/comm_mpi_particle.hpp"
|
---|
[e9cfc4] | 73 |
|
---|
[0990e1] | 74 | #include "CodePatterns/MemDebug.hpp"
|
---|
| 75 |
|
---|
| 76 | #include "Jobs/VMGJob.hpp"
|
---|
| 77 |
|
---|
[e9cfc4] | 78 | #include "LinearAlgebra/defs.hpp"
|
---|
| 79 | #include "Jobs/InterfaceVMGJob.hpp"
|
---|
| 80 |
|
---|
[a82602] | 81 | #include "CodePatterns/Assert.hpp"
|
---|
| 82 |
|
---|
[e9cfc4] | 83 | using namespace VMG;
|
---|
[d2a0f6d] | 84 |
|
---|
[d12d621] | 85 | VMGJob::VMGJob(
|
---|
| 86 | const JobId_t _JobId,
|
---|
[cd77fc] | 87 | const SamplingGrid &_density_grid,
|
---|
| 88 | const std::vector< std::vector< double > > &_particle_positions,
|
---|
[065574] | 89 | const std::vector< double > &_particle_charges,
|
---|
[cd2591] | 90 | const size_t _near_field_cells,
|
---|
[b6b21a] | 91 | const size_t _interpolation_degree,
|
---|
[e2925fd] | 92 | const bool _DoImportParticles,
|
---|
[ee9018] | 93 | const bool _DoPrintDebug,
|
---|
| 94 | const bool _OpenBoundaryConditions) :
|
---|
[d2a0f6d] | 95 | FragmentJob(_JobId),
|
---|
[442cee] | 96 | density_grid(_density_grid),
|
---|
[cd77fc] | 97 | particle_positions(_particle_positions),
|
---|
[2bc560] | 98 | particle_charges(_particle_charges),
|
---|
[065574] | 99 | near_field_cells(_near_field_cells),
|
---|
[cd2591] | 100 | interpolation_degree(_interpolation_degree),
|
---|
[e2925fd] | 101 | DoImportParticles(_DoImportParticles),
|
---|
| 102 | DoPrintDebug(_DoPrintDebug),
|
---|
[ee9018] | 103 | OpenBoundaryConditions(_OpenBoundaryConditions),
|
---|
[a82602] | 104 | returndata(static_cast<const SamplingGridProperties &>(_density_grid)),
|
---|
[e2925fd] | 105 | particles()
|
---|
[28c025] | 106 | {}
|
---|
[d2a0f6d] | 107 |
|
---|
| 108 | VMGJob::VMGJob() :
|
---|
[065574] | 109 | FragmentJob(JobId::IllegalJob),
|
---|
[cd2591] | 110 | near_field_cells(3),
|
---|
| 111 | interpolation_degree(3),
|
---|
[e2925fd] | 112 | DoImportParticles(true),
|
---|
| 113 | DoPrintDebug(false),
|
---|
[ee9018] | 114 | OpenBoundaryConditions(false),
|
---|
[e2925fd] | 115 | particles()
|
---|
[e089fb] | 116 | {}
|
---|
| 117 |
|
---|
| 118 | VMGJob::~VMGJob()
|
---|
| 119 | {
|
---|
| 120 | }
|
---|
| 121 |
|
---|
| 122 |
|
---|
| 123 | VMGJob::particle_arrays::particle_arrays() :
|
---|
[a82602] | 124 | num_particles(0),
|
---|
| 125 | f(NULL),
|
---|
| 126 | x(NULL),
|
---|
| 127 | p(NULL),
|
---|
| 128 | q(NULL)
|
---|
[d2a0f6d] | 129 | {}
|
---|
| 130 |
|
---|
[e089fb] | 131 | VMGJob::particle_arrays::~particle_arrays()
|
---|
[a82602] | 132 | {
|
---|
| 133 | delete[] f;
|
---|
| 134 | delete[] x;
|
---|
| 135 | delete[] p;
|
---|
| 136 | delete[] q;
|
---|
| 137 | }
|
---|
| 138 |
|
---|
[e089fb] | 139 | void VMGJob::particle_arrays::init(const size_t _num_particles)
|
---|
[a82602] | 140 | {
|
---|
[e089fb] | 141 | num_particles = _num_particles;
|
---|
[a82602] | 142 | f = new double[num_particles*3];
|
---|
| 143 | x = new double[num_particles*3];
|
---|
| 144 | p = new double[num_particles];
|
---|
| 145 | q = new double[num_particles];
|
---|
[e089fb] | 146 | }
|
---|
| 147 |
|
---|
| 148 | void VMGJob::InitVMGArrays()
|
---|
| 149 | {
|
---|
| 150 | particles.init(particle_charges.size());
|
---|
[a82602] | 151 |
|
---|
| 152 | {
|
---|
| 153 | size_t index=0;
|
---|
| 154 | for (std::vector< std::vector< double > >::const_iterator iter = particle_positions.begin();
|
---|
| 155 | iter != particle_positions.end(); ++iter) {
|
---|
| 156 | for (std::vector< double >::const_iterator positer = (*iter).begin();
|
---|
| 157 | positer != (*iter).end(); ++positer) {
|
---|
[e089fb] | 158 | particles.f[index] = 0.;
|
---|
| 159 | particles.x[index++] = *positer;
|
---|
[a82602] | 160 | }
|
---|
| 161 | }
|
---|
[e089fb] | 162 | ASSERT( index == particles.num_particles*3,
|
---|
[a82602] | 163 | "VMGJob::VMGJob() - too many particles or components in particle_positions.");
|
---|
| 164 | }
|
---|
| 165 | {
|
---|
| 166 | size_t index = 0;
|
---|
| 167 | for (std::vector< double >::const_iterator iter = particle_charges.begin();
|
---|
| 168 | iter != particle_charges.end(); ++iter) {
|
---|
[e089fb] | 169 | particles.p[index] = 0.;
|
---|
| 170 | particles.q[index++] = *iter;
|
---|
[a82602] | 171 | }
|
---|
[e089fb] | 172 | ASSERT( index == particles.num_particles,
|
---|
[a82602] | 173 | "VMGJob::VMGJob() - too many charges in particle_charges.");
|
---|
| 174 | }
|
---|
| 175 | }
|
---|
| 176 |
|
---|
[d2a0f6d] | 177 |
|
---|
| 178 | FragmentResult::ptr VMGJob::Work()
|
---|
| 179 | {
|
---|
[e9cfc4] | 180 | // initialize VMG library solver
|
---|
| 181 | InitVMG();
|
---|
| 182 |
|
---|
| 183 | /*
|
---|
| 184 | * Start the multigrid solver
|
---|
| 185 | */
|
---|
| 186 | MG::Solve();
|
---|
| 187 |
|
---|
| 188 | /// create and fill result object
|
---|
[442cee] | 189 | // place into returnstream
|
---|
| 190 | std::stringstream returnstream;
|
---|
[ae68b5] | 191 | {
|
---|
| 192 | const VMGData archivedata(returndata);
|
---|
| 193 | boost::archive::text_oarchive oa(returnstream);
|
---|
| 194 | oa << archivedata;
|
---|
| 195 | }
|
---|
[442cee] | 196 |
|
---|
| 197 | FragmentResult::ptr ptr( new FragmentResult(getId(), returnstream.str()) );
|
---|
[e9cfc4] | 198 |
|
---|
| 199 | /*
|
---|
| 200 | * Delete all data.
|
---|
| 201 | */
|
---|
| 202 | MG::Destroy();
|
---|
| 203 |
|
---|
[d2a0f6d] | 204 | return ptr;
|
---|
| 205 | }
|
---|
| 206 |
|
---|
[e9cfc4] | 207 | /** Initialization of VMG library.
|
---|
| 208 | *
|
---|
| 209 | * The contents is heavily inspired from interface_fcs.cpp: VMG_fcs_init() of
|
---|
| 210 | * the ScaFaCoS project.
|
---|
| 211 | *
|
---|
| 212 | */
|
---|
| 213 | void VMGJob::InitVMG()
|
---|
| 214 | {
|
---|
[ee9018] | 215 | Boundary *boundary = NULL;
|
---|
| 216 | if (OpenBoundaryConditions)
|
---|
| 217 | boundary = new Boundary(Open, Open, Open);
|
---|
| 218 | else
|
---|
| 219 | boundary = new Boundary(Periodic, Periodic, Periodic);
|
---|
[e9cfc4] | 220 |
|
---|
| 221 | /*
|
---|
| 222 | * Choose multigrid components (self-registering)
|
---|
| 223 | */
|
---|
[ee9018] | 224 |
|
---|
[8a8c8c] | 225 | #ifdef HAVE_MPI
|
---|
[ee9018] | 226 | new Particle::CommMPI(*boundary, new DomainDecompositionMPI());
|
---|
[8a8c8c] | 227 | #else
|
---|
[ee9018] | 228 | new CommSerial(*boundary);
|
---|
[8a8c8c] | 229 | #endif
|
---|
[ee9018] | 230 | if (OpenBoundaryConditions)
|
---|
| 231 | new DiscretizationPoissonFV(2);
|
---|
| 232 | else
|
---|
| 233 | new DiscretizationPoissonFD(4);
|
---|
[e9cfc4] | 234 | new VMGInterfaces::InterfaceVMGJob(
|
---|
[8f3cdd] | 235 | density_grid,
|
---|
[2bc560] | 236 | returndata,
|
---|
[cd77fc] | 237 | particle_positions,
|
---|
| 238 | particle_charges,
|
---|
[ee9018] | 239 | *boundary,
|
---|
[e9cfc4] | 240 | 2,
|
---|
[28c025] | 241 | density_grid.level,
|
---|
| 242 | Vector(density_grid.begin),
|
---|
[3d9a8d] | 243 | density_grid.end[0]-density_grid.begin[0],
|
---|
[b6b21a] | 244 | near_field_cells,
|
---|
[e2925fd] | 245 | DoImportParticles ?
|
---|
| 246 | VMGInterfaces::InterfaceVMGJob::DoImportParticles
|
---|
| 247 | : VMGInterfaces::InterfaceVMGJob::DontImportParticles,
|
---|
[b6b21a] | 248 | DoPrintDebug);
|
---|
[e9cfc4] | 249 | const int cycle_type = 1; // V-type
|
---|
[ee9018] | 250 | if (OpenBoundaryConditions) {
|
---|
| 251 | new LevelOperatorFAS(Stencils::RestrictionFullWeight, Stencils::Injection, Stencils::InterpolationTrilinear);
|
---|
| 252 | new Givens<SolverRegular>();
|
---|
| 253 | new CycleFASDirichlet(cycle_type);
|
---|
| 254 | new GaussSeidelRBPoisson2();
|
---|
| 255 | } else {
|
---|
| 256 | new LevelOperatorCS(Stencils::RestrictionFullWeight, Stencils::InterpolationTrilinear);
|
---|
| 257 | new Givens<SolverSingular>();
|
---|
| 258 | new CycleCSPeriodic(cycle_type);
|
---|
| 259 | new GaussSeidelRBPoisson4();
|
---|
| 260 | }
|
---|
| 261 | delete boundary;
|
---|
[e9cfc4] | 262 |
|
---|
| 263 | /*
|
---|
| 264 | * Register required parameters
|
---|
| 265 | */
|
---|
| 266 | new ObjectStorage<int>("PRESMOOTHSTEPS", 3);
|
---|
| 267 | new ObjectStorage<int>("POSTSMOOTHSTEPS", 3);
|
---|
| 268 | new ObjectStorage<vmg_float>("PRECISION", 1.0e-10);
|
---|
| 269 | new ObjectStorage<int>("MAX_ITERATION", 15);
|
---|
[cd77fc] | 270 | new ObjectStorage<int>("PARTICLE_NEAR_FIELD_CELLS", near_field_cells);
|
---|
[cd2591] | 271 | new ObjectStorage<int>("PARTICLE_INTERPOLATION_DEGREE", interpolation_degree);
|
---|
[e9cfc4] | 272 |
|
---|
[a82602] | 273 | // first initialize f,x,p,q,... from STL vectors
|
---|
| 274 | InitVMGArrays();
|
---|
| 275 | // then initialize as objects with VMG's storage
|
---|
[e089fb] | 276 | new ObjectStorage<vmg_float*>("PARTICLE_POS_ARRAY", particles.x);
|
---|
| 277 | new ObjectStorage<vmg_float*>("PARTICLE_CHARGE_ARRAY", particles.q);
|
---|
| 278 | new ObjectStorage<vmg_float*>("PARTICLE_POTENTIAL_ARRAY", particles.p);
|
---|
| 279 | new ObjectStorage<vmg_float*>("PARTICLE_FIELD_ARRAY", particles.f);
|
---|
| 280 | new ObjectStorage<vmg_int>("PARTICLE_NUM_LOCAL", particles.num_particles);
|
---|
[a82602] | 281 |
|
---|
[e9cfc4] | 282 | /*
|
---|
| 283 | * Post init
|
---|
| 284 | */
|
---|
| 285 | MG::PostInit();
|
---|
| 286 |
|
---|
| 287 | /*
|
---|
| 288 | * Check whether the library is correctly initialized now.
|
---|
| 289 | */
|
---|
| 290 | MG::IsInitialized();
|
---|
| 291 | }
|
---|
| 292 |
|
---|
[d2a0f6d] | 293 | // we need to explicitly instantiate the serialization functions as
|
---|
| 294 | // its is only serialized through its base class FragmentJob
|
---|
| 295 | BOOST_CLASS_EXPORT_IMPLEMENT(VMGJob)
|
---|