DFT-FE 1.1.0-pre
Density Functional Theory With Finite-Elements
Loading...
Searching...
No Matches
force.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (c) 2017-2025 The Regents of the University of Michigan and DFT-FE
4// authors.
5//
6// This file is part of the DFT-FE code.
7//
8// The DFT-FE code is free software; you can use it, redistribute
9// it, and/or modify it under the terms of the GNU Lesser General
10// Public License as published by the Free Software Foundation; either
11// version 2.1 of the License, or (at your option) any later version.
12// The full text of the license can be found in the file LICENSE at
13// the top level of the DFT-FE distribution.
14//
15// ---------------------------------------------------------------------
16
17
18#ifndef force_H_
19#define force_H_
20#include "vselfBinsManager.h"
21#include "dftParameters.h"
22
23#include "constants.h"
24#include "headers.h"
26#include <dftd.h>
27#include <oncvClass.h>
29#include <feevaluationWrapper.h>
30
31namespace dftfe
32{
33 // forward declaration
34 template <dftfe::utils::MemorySpace memory>
35 class dftClass;
36
37 /**
38 * @brief computes configurational forces in KSDFT
39 *
40 * This class computes and stores the configurational forces corresponding to
41 * geometry optimization. It uses the formulation in the paper by Motamarri
42 * et.al. (https://link.aps.org/doi/10.1103/PhysRevB.97.165132) which provides
43 * an unified approach to atomic forces corresponding to internal atomic
44 * relaxation and cell stress corresponding to cell relaxation.
45 *
46 * @author Sambit Das
47 */
48 template <dftfe::utils::MemorySpace memorySpace>
50 {
51 friend class dftClass<memorySpace>;
52
53 public:
54 /** @brief Constructor.
55 *
56 * @param _dftPtr pointer to dftClass
57 * @param mpi_comm_parent parent mpi_communicator
58 * @param mpi_comm_domain domain decomposition mpi_communicator
59 */
61 const MPI_Comm &mpi_comm_parent,
62 const MPI_Comm &mpi_comm_domain,
63 const dftParameters &dftParams);
64
65 /** @brief initializes data structures inside forceClass assuming unmoved triangulation.
66 *
67 * initUnmoved is the first step of the initialization/reinitialization of
68 * force class when starting from a new unmoved triangulation. It creates
69 * the dofHandler with linear finite elements and three components
70 * corresponding to the three force components. It also creates the
71 * corresponding constraint matrices which is why an unmoved triangulation
72 * is necessary. Finally this function also initializes the gaussianMovePar
73 * data member.
74 *
75 * @param triangulation reference to unmoved triangulation where the mesh nodes have not
76 * been manually moved.
77 * @param isElectrostaticsMesh boolean parameter specifying whether this triangulatio is to be used for
78 * for the electrostatics part of the configurational force.
79 * @return void.
80 */
81 void
82 initUnmoved(const dealii::Triangulation<3, 3> &triangulation,
83 const dealii::Triangulation<3, 3> &serialTriangulation,
84 const std::vector<std::vector<double>> &domainBoundingVectors,
85 const bool isElectrostaticsMesh);
86
87 /** @brief initializes data structures inside forceClass which depend on the moved mesh.
88 *
89 * initMoved is the second step (first step call initUnmoved) of the
90 * initialization/reinitialization of force class when starting from a new
91 * mesh, and the first step when recomputing forces on a perturbed mesh.
92 * initMoved assumes that the triangulation whose reference was passed to
93 * the forceClass object in the initUnmoved call now has its nodes moved
94 * such that all atomic positions lie on nodes.
95 *
96 * @return void.
97 */
98 void
100 std::vector<const dealii::DoFHandler<3> *> &dofHandlerVectorMatrixFree,
101 std::vector<const dealii::AffineConstraints<double> *>
102 &constraintsVectorMatrixFree,
103 const bool isElectrostaticsMesh);
104
105 /** @brief initializes and precomputes pseudopotential related data structuers required for configurational force
106 * and stress computation.
107 *
108 * This function is only activated for pseudopotential calculations and is
109 * currently called when initializing/reinitializing the dftClass object.
110 * This function initializes and precomputes the pseudopotential
111 * datastructures for local and non-local parts. Separate internal function
112 * calls are made for KB and ONCV projectors.
113 *
114 * @return void.
115 */
116 void
118
119 /** @brief computes the configurational force on all atoms corresponding to a Gaussian generator,
120 * which represents perturbation of the underlying space.
121 *
122 * The Gaussian generator is taken to be exp(-d_gaussianConstant*r^2), r
123 * being the distance from the atom. Currently d_gaussianConstant is
124 * hardcoded to be 4.0. To get the computed atomic forces use getAtomsForces
125 *
126 * @return void.
127 */
128 void
130 const dealii::MatrixFree<3, double> &matrixFreeData,
131 const dispersionCorrection &dispersionCorr,
132 const dftfe::uInt eigenDofHandlerIndex,
133 const dftfe::uInt smearedChargeQuadratureId,
134 const dftfe::uInt lpspQuadratureIdElectro,
135 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
136 const dftfe::uInt phiTotDofHandlerIndexElectro,
137 const distributedCPUVec<double> &phiTotRhoOutElectro,
138 const std::vector<
140 &rhoOutValues,
141 const std::vector<
143 &gradRhoOutValues,
145 &rhoTotalOutValuesLpsp,
147 &gradRhoTotalOutValuesLpsp,
148 const std::map<dealii::CellId, std::vector<double>> &rhoCoreValues,
149 const std::map<dealii::CellId, std::vector<double>> &gradRhoCoreValues,
150 const std::map<dealii::CellId, std::vector<double>> &hessianRhoCoreValues,
151 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
152 &gradRhoCoreAtoms,
153 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
154 &hessianRhoCoreAtoms,
155 const std::map<dealii::CellId, std::vector<double>> &pseudoVLocElectro,
156 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
157 &pseudoVLocAtomsElectro,
158 const dealii::AffineConstraints<double> &hangingPlusPBCConstraintsElectro,
159 const vselfBinsManager &vselfBinsManagerElectro);
160
161 /** @brief returns a copy of the configurational force on all global atoms.
162 *
163 * computeAtomsForces must be called prior to this function call.
164 *
165 * @return std::vector<double> flattened array of the configurational force on all atoms,
166 * the three force components on each atom being the leading dimension.
167 * Units- Hartree/Bohr
168 */
169 std::vector<double> &
171
172 /** @brief prints the currently stored configurational forces on atoms and the Gaussian generator constant
173 * used to compute them.
174 *
175 * @return void.
176 */
177 void
179
180 /** @brief Update force generator Gaussian constant.
181 *
182 * @return void.
183 */
184 // void updateGaussianConstant(const double newGaussianConstant);
185
186 /** @brief computes the configurational stress on the domain corresponding to
187 * affine deformation of the periodic cell.
188 *
189 * This function cannot be called for fully non-periodic calculations.
190 *
191 * @return void.
192 */
193 void
195 const dealii::MatrixFree<3, double> &matrixFreeData,
196 const dispersionCorrection &dispersionCorr,
197 const dftfe::uInt eigenDofHandlerIndex,
198 const dftfe::uInt smearedChargeQuadratureId,
199 const dftfe::uInt lpspQuadratureIdElectro,
200 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
201 const dftfe::uInt phiTotDofHandlerIndexElectro,
202 const distributedCPUVec<double> &phiTotRhoOutElectro,
203 const std::vector<
205 &rhoOutValues,
206 const std::vector<
208 &gradRhoOutValues,
210 &rhoTotalOutValuesLpsp,
212 &gradRhoTotalOutValuesLpsp,
213 const std::map<dealii::CellId, std::vector<double>> &pseudoVLocElectro,
214 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
215 &pseudoVLocAtomsElectro,
216 const std::map<dealii::CellId, std::vector<double>> &rhoCoreValues,
217 const std::map<dealii::CellId, std::vector<double>> &gradRhoCoreValues,
218 const std::map<dealii::CellId, std::vector<double>> &hessianRhoCoreValues,
219 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
220 &gradRhoCoreAtoms,
221 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
222 &hessianRhoCoreAtoms,
223 const dealii::AffineConstraints<double> &hangingPlusPBCConstraintsElectro,
224 const vselfBinsManager &vselfBinsManagerElectro);
225
226 /** @brief prints the currently stored configurational stress tensor.
227 *
228 * @return void.
229 */
230 void
232
233 /** @brief returns a copy of the current stress tensor value.
234 *
235 * computeStress must be call prior to this function call.
236 *
237 * @return dealii::Tensor<2,3,double> second order stress Tensor in Hartree/Bohr^3
238 */
239 dealii::Tensor<2, 3, double> &
241
242 /** @brief get the value of Gaussian generator parameter (d_gaussianConstant).
243 * Gaussian generator: Gamma(r)= exp(-d_gaussianConstant*r^2).
244 *
245 */
246 // double getGaussianGeneratorParameter() const;
247
248 private:
249 /** @brief Locates and stores the global dof indices of d_dofHandlerForce whose cooridinates match
250 * with the atomic positions.
251 *
252 * @return void.
253 */
254 void
255 locateAtomCoreNodesForce(const dealii::DoFHandler<3> &dofHandlerForce,
256 const dealii::IndexSet &locally_owned_dofsForce,
257 std::map<std::pair<dftfe::uInt, dftfe::uInt>,
258 dftfe::uInt> &atomsForceDofs);
259
260 void
262 const dealii::DoFHandler<3> &dofHandler,
263 const dealii::DoFHandler<3> &dofHandlerForce,
264 const dealii::AffineConstraints<double> &hangingPlusPBCConstraints,
266 std::vector<std::vector<dealii::DoFHandler<3>::active_cell_iterator>>
267 &cellsVselfBallsDofHandler,
268 std::vector<std::vector<dealii::DoFHandler<3>::active_cell_iterator>>
269 &cellsVselfBallsDofHandlerForce,
270 std::vector<std::map<dealii::CellId, dftfe::uInt>>
271 &cellsVselfBallsClosestAtomIdDofHandler,
272 std::map<dftfe::uInt, dftfe::uInt> &AtomIdBinIdLocalDofHandler,
273 std::vector<std::map<dealii::DoFHandler<3>::active_cell_iterator,
274 std::vector<dftfe::uInt>>>
275 &cellFacesVselfBallSurfacesDofHandler,
276 std::vector<std::map<dealii::DoFHandler<3>::active_cell_iterator,
277 std::vector<dftfe::uInt>>>
278 &cellFacesVselfBallSurfacesDofHandlerForce);
279
280 void
282 const dealii::MatrixFree<3, double> &matrixFreeData,
283 const dealii::MatrixFree<3, double> &matrixFreeDataElectro);
284
285 void
287
288 void
290 const dealii::MatrixFree<3, double> &matrixFreeData,
291 const dftfe::uInt eigenDofHandlerIndex,
292 const dftfe::uInt smearedChargeQuadratureId,
293 const dftfe::uInt lpspQuadratureIdElectro,
294 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
295 const dftfe::uInt phiTotDofHandlerIndexElectro,
296 const distributedCPUVec<double> &phiTotRhoOutElectro,
297 const std::vector<
299 &rhoOutValues,
300 const std::vector<
302 &gradRhoOutValues,
304 &rhoTotalOutValuesLpsp,
306 &gradRhoTotalOutValuesLpsp,
307 const std::map<dealii::CellId, std::vector<double>> &rhoCoreValues,
308 const std::map<dealii::CellId, std::vector<double>> &gradRhoCoreValues,
309 const std::map<dealii::CellId, std::vector<double>> &hessianRhoCoreValues,
310 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
311 &gradRhoCoreAtoms,
312 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
313 &hessianRhoCoreAtoms,
314 const std::map<dealii::CellId, std::vector<double>> &pseudoVLocElectro,
315 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
316 &pseudoVLocAtomsElectro,
317 const vselfBinsManager &vselfBinsManagerElectro);
318
319 void
321 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
322 const dftfe::uInt phiTotDofHandlerIndexElectro,
323 const dftfe::uInt smearedChargeQuadratureId,
324 const dftfe::uInt lpspQuadratureIdElectro,
325 const distributedCPUVec<double> &phiTotRhoOutElectro,
327 &rhoTotalOutValues,
329 &rhoTotalOutValuesLpsp,
331 &gradRhoTotalOutValues,
333 &gradRhoTotalOutValuesLpsp,
334 const std::map<dealii::CellId, std::vector<double>> &pseudoVLocElectro,
335 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
336 &pseudoVLocAtomsElectro,
337 const vselfBinsManager &vselfBinsManagerElectro);
338
339 void
341
342 void
344 const dealii::DoFHandler<3> &dofHandlerElectro,
345 const vselfBinsManager &vselfBinsManagerElectro,
346 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
347 const dftfe::uInt smearedChargeQuadratureId);
348
349 void
351
352 void
354 const dealii::MatrixFree<3, double> &matrixFreeData,
355 const dftfe::uInt eigenDofHandlerIndex,
356 const dftfe::uInt smearedChargeQuadratureId,
357 const dftfe::uInt lpspQuadratureIdElectro,
358 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
359 const dftfe::uInt phiTotDofHandlerIndexElectro,
360 const distributedCPUVec<double> &phiTotRhoOutElectro,
361 const std::vector<
363 &rhoOutValues,
364 const std::vector<
366 &gradRhoOutValues,
368 &rhoTotalOutValuesLpsp,
370 &gradRhoTotalOutValuesLpsp,
371 const std::map<dealii::CellId, std::vector<double>> &rhoCoreValues,
372 const std::map<dealii::CellId, std::vector<double>> &gradRhoCoreValues,
373 const std::map<dealii::CellId, std::vector<double>> &hessianRhoCoreValues,
374 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
375 &gradRhoCoreAtoms,
376 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
377 &hessianRhoCoreAtoms,
378 const std::map<dealii::CellId, std::vector<double>> &pseudoVLocElectro,
379 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
380 &pseudoVLocAtomsElectro,
381 const vselfBinsManager &vselfBinsManagerElectro);
382
383 void
385 std::map<dftfe::uInt, std::vector<double>>
386 &forceContributionFPSPLocalGammaAtoms,
387 dealii::FEValues<3> &feValues,
388 dealii::FEFaceValues<3> &feFaceValues,
390 const dealii::MatrixFree<3, double> &matrixFreeData,
391 const dftfe::uInt phiTotDofHandlerIndexElectro,
392 const dftfe::uInt cell,
393 const dealii::AlignedVector<dealii::VectorizedArray<double>> &rhoQuads,
394 const dealii::AlignedVector<
395 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>> &gradRhoQuads,
396 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
397 &pseudoVLocAtoms,
399 const std::vector<std::map<dealii::CellId, dftfe::uInt>>
400 &cellsVselfBallsClosestAtomIdDofHandler);
401
402 void
404 std::map<dftfe::uInt, std::vector<double>>
405 &forceContributionSmearedChargesGammaAtoms,
407 const dealii::MatrixFree<3, double> &matrixFreeData,
408 const dftfe::uInt cell,
409 const dealii::AlignedVector<
410 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>> &gradPhiTotQuads,
411 const std::vector<dftfe::uInt> &nonTrivialAtomIdsMacroCell,
412 const std::map<dealii::CellId, std::vector<dftfe::Int>>
413 &bQuadAtomIdsAllAtoms,
414 const dealii::AlignedVector<dealii::VectorizedArray<double>>
415 &smearedbQuads);
416
417 void
419 std::map<dftfe::uInt, std::vector<double>>
420 &forceContributionSmearedChargesGammaAtoms,
422 const dealii::MatrixFree<3, double> &matrixFreeData,
423 const dftfe::uInt cell,
424 const dealii::AlignedVector<
425 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>>
426 &gradVselfBinQuads,
427 const std::vector<dftfe::uInt> &nonTrivialAtomIdsMacroCell,
428 const std::map<dealii::CellId, std::vector<dftfe::Int>>
429 &bQuadAtomIdsAllAtoms,
430 const dealii::AlignedVector<dealii::VectorizedArray<double>>
431 &smearedbQuads);
432
433 void
435 std::map<dftfe::uInt, std::vector<double>>
436 &forceContributionFNonlinearCoreCorrectionGammaAtoms,
438 const dealii::MatrixFree<3, double> &matrixFreeData,
439 const dftfe::uInt cell,
440 const dealii::AlignedVector<dealii::VectorizedArray<double>> &vxcQuads,
441 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
442 &gradRhoCoreAtoms);
443
444
445 void
447 std::map<dftfe::uInt, std::vector<double>>
448 &forceContributionFNonlinearCoreCorrectionGammaAtoms,
450 const dealii::MatrixFree<3, double> &matrixFreeData,
451 const dftfe::uInt cell,
452 const dealii::AlignedVector<
453 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>> &derExcGradRho,
454 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
455 &hessianRhoCoreAtoms);
456
457 void
459 std::map<dftfe::uInt, std::vector<double>>
460 &forceContributionFNonlinearCoreCorrectionGammaAtoms,
462 const dealii::MatrixFree<3, double> &matrixFreeData,
463 const dftfe::uInt cell,
464 const dealii::AlignedVector<dealii::VectorizedArray<double>>
465 &vxcQuadsSpin0,
466 const dealii::AlignedVector<dealii::VectorizedArray<double>>
467 &vxcQuadsSpin1,
468 const dealii::AlignedVector<
469 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>>
470 &derExcGradRhoSpin0,
471 const dealii::AlignedVector<
472 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>>
473 &derExcGradRhoSpin1,
474 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
475 &gradRhoCoreAtoms,
476 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
477 &hessianRhoCoreAtoms,
478 const bool isXCGGA = false);
479
480 void
482 const std::map<dftfe::uInt, std::vector<double>>
483 &forceContributionFPSPLocalGammaAtoms,
484 const std::map<std::pair<dftfe::uInt, dftfe::uInt>, dftfe::uInt>
485 &atomsForceDofs,
486 const dealii::AffineConstraints<double> &constraintsNoneForce,
487 distributedCPUVec<double> &configForceVectorLinFE);
488
489 void
491 const std::map<dftfe::uInt, std::vector<double>>
492 &forceContributionLocalGammaAtoms,
493 std::vector<double> &accumForcesVector);
494
495
496 void
498 std::map<dftfe::uInt, std::vector<double>>
499 &forceContributionFnlGammaAtoms,
500 const dealii::MatrixFree<3, double> &matrixFreeData,
501 FEEvaluationWrapperClass<3> &forceEvalNLP,
502 const std::shared_ptr<
504 nonLocalOp,
505 dftfe::uInt numNonLocalAtomsCurrentProcess,
506 const std::vector<dftfe::Int> &globalChargeIdNonLocalAtoms,
507 const std::vector<dftfe::uInt> &numberPseudoWaveFunctionsPerAtom,
508 const dftfe::uInt cell,
509 const std::map<dealii::CellId, dftfe::uInt> &cellIdToCellNumberMap,
510#ifdef USE_COMPLEX
511 const std::vector<dataTypes::number>
512 &projectorKetTimesPsiTimesVTimesPartOccContractionPsiQuadsFlattened,
513#endif
514 const std::vector<dataTypes::number> &zetaDeltaVQuadsFlattened,
515 const std::vector<dataTypes::number> &
516 projectorKetTimesPsiTimesVTimesPartOccContractionGradPsiQuadsFlattened);
517
518
519 void
521 dealii::AlignedVector<
522 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>> &FVectQuads,
523 const dealii::MatrixFree<3, double> &matrixFreeData,
524 const dftfe::uInt numQuadPoints,
525 const std::shared_ptr<
527 nonLocalOp,
528 const dftfe::uInt numNonLocalAtomsCurrentProcess,
529 const std::vector<dftfe::Int> &globalChargeIdNonLocalAtoms,
530 const std::vector<dftfe::uInt> &numberPseudoWaveFunctionsPerAtom,
531 const dftfe::uInt cell,
532 const std::map<dealii::CellId, dftfe::uInt> &cellIdToCellNumberMap,
533 const std::vector<dataTypes::number> &zetaDeltaVQuadsFlattened,
534 const std::vector<dataTypes::number> &
535 projectorKetTimesPsiTimesVTimesPartOccContractionGradPsiQuadsFlattened);
536
537 void
539 const std::map<dftfe::uInt, std::vector<double>>
540 &forceContributionFnlGammaAtoms);
541
542 void
544 dealii::Tensor<2, 3, double> &stressContribution,
545 const dealii::MatrixFree<3, double> &matrixFreeData,
546 const dftfe::uInt numQuadPoints,
547 const std::vector<double> &jxwQuadsSubCells,
548 const dftfe::uInt cell,
549 const dftfe::uInt numNonLocalAtomsCurrentProcess,
550 const std::shared_ptr<
552 nonLocalOp,
553 const std::vector<dftfe::uInt> &numberPseudoWaveFunctionsPerAtom,
554 const std::map<dealii::CellId, dftfe::uInt> &cellIdToCellNumberMap,
555 const std::vector<dataTypes::number> &zetalmDeltaVlProductDistImageAtoms,
556#ifdef USE_COMPLEX
557 const std::vector<dataTypes::number>
558 &projectorKetTimesPsiTimesVTimesPartOccContractionPsiQuadsFlattened,
559#endif
560 const std::vector<dataTypes::number>
561 &projectorKetTimesPsiTimesVTimesPartOccContractionGradPsiQuadsFlattened,
562 const bool isSpinPolarized);
563
564 void
566 bool allowGaussianOverlapOnAtoms = false);
567
568 void
570
571 void
573 const dealii::DoFHandler<3> &dofHandlerElectro,
574 const vselfBinsManager &vselfBinsManagerElectro,
575 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
576 const dftfe::uInt smearedChargeQuadratureId);
577
578 void
580 const dealii::MatrixFree<3, double> &matrixFreeData,
581 const dftfe::uInt eigenDofHandlerIndex,
582 const dftfe::uInt smearedChargeQuadratureId,
583 const dftfe::uInt lpspQuadratureIdElectro,
584 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
585 const dftfe::uInt phiTotDofHandlerIndexElectro,
586 const distributedCPUVec<double> &phiTotRhoOutElectro,
587 const std::vector<
589 &rhoOutValues,
590 const std::vector<
592 &gradRhoOutValues,
594 &rhoTotalOutValuesLpsp,
596 &gradRhoTotalOutValuesLpsp,
597 const std::map<dealii::CellId, std::vector<double>> &pseudoVLocElectro,
598 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
599 &pseudoVLocAtomsElectro,
600 const std::map<dealii::CellId, std::vector<double>> &rhoCoreValues,
601 const std::map<dealii::CellId, std::vector<double>> &gradRhoCoreValues,
602 const std::map<dealii::CellId, std::vector<double>> &hessianRhoCoreValues,
603 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
604 &gradRhoCoreAtoms,
605 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
606 &hessianRhoCoreAtoms,
607 const vselfBinsManager &vselfBinsManagerElectro);
608
609 void
611 const dealii::MatrixFree<3, double> &matrixFreeDataElectro,
612 const dftfe::uInt phiTotDofHandlerIndexElectro,
613 const dftfe::uInt smearedChargeQuadratureId,
614 const dftfe::uInt lpspQuadratureIdElectro,
615 const distributedCPUVec<double> &phiTotRhoOutElectro,
617 &rhoTotalOutValues,
619 &rhoTotalOutValuesLpsp,
621 &gradRhoTotalOutValuesElectro,
623 &gradRhoTotalOutValuesElectroLpsp,
624 const std::map<dealii::CellId, std::vector<double>> &pseudoVLocElectro,
625 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
626 &pseudoVLocAtomsElectro,
627 const vselfBinsManager &vselfBinsManagerElectro);
628
629 void
631 dealii::FEValues<3> &feValues,
632 dealii::FEFaceValues<3> &feFaceValues,
634 const dealii::MatrixFree<3, double> &matrixFreeData,
635 const dftfe::uInt phiTotDofHandlerIndexElectro,
636 const dftfe::uInt cell,
637 const dealii::AlignedVector<dealii::VectorizedArray<double>> &rhoQuads,
638 const dealii::AlignedVector<
639 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>> &gradRhoQuads,
640 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
641 &pseudoVLocAtoms,
643 const std::vector<std::map<dealii::CellId, dftfe::uInt>>
644 &cellsVselfBallsClosestAtomIdDofHandler);
645
646 void
649 const dealii::MatrixFree<3, double> &matrixFreeData,
650 const dftfe::uInt cell,
651 const dealii::AlignedVector<dealii::VectorizedArray<double>> &vxcQuads,
652 const dealii::AlignedVector<
653 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>> &derExcGradRho,
654 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
655 &gradRhoCoreAtoms,
656 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
657 &hessianRhoCoreAtoms);
658
659 void
662 const dealii::MatrixFree<3, double> &matrixFreeData,
663 const dftfe::uInt cell,
664 const dealii::AlignedVector<dealii::VectorizedArray<double>>
665 &vxcQuadsSpin0,
666 const dealii::AlignedVector<dealii::VectorizedArray<double>>
667 &vxcQuadsSpin1,
668 const dealii::AlignedVector<
669 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>>
670 &derExcGradRhoSpin0,
671 const dealii::AlignedVector<
672 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>>
673 &derExcGradRhoSpin1,
674 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
675 &gradRhoCoreAtoms,
676 const std::map<dftfe::uInt, std::map<dealii::CellId, std::vector<double>>>
677 &hessianRhoCoreAtoms,
678 const bool isXCGGA = false);
679
680 void
683 const dealii::MatrixFree<3, double> &matrixFreeData,
684 const dftfe::uInt cell,
685 const dealii::AlignedVector<
686 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>> &gradPhiTotQuads,
687 const std::vector<dftfe::uInt> &nonTrivialAtomImageIdsMacroCell,
688 const std::map<dealii::CellId, std::vector<dftfe::Int>>
689 &bQuadAtomIdsAllAtomsImages,
690 const dealii::AlignedVector<dealii::VectorizedArray<double>>
691 &smearedbQuads);
692
693 void
696 const dealii::MatrixFree<3, double> &matrixFreeData,
697 const dftfe::uInt cell,
698 const dealii::AlignedVector<
699 dealii::Tensor<1, 3, dealii::VectorizedArray<double>>> &gradVselfQuads,
700 const std::vector<dftfe::uInt> &nonTrivialAtomImageIdsMacroCell,
701 const std::map<dealii::CellId, std::vector<dftfe::Int>>
702 &bQuadAtomIdsAllAtomsImages,
703 const dealii::AlignedVector<dealii::VectorizedArray<double>>
704 &smearedbQuads);
705
706 void
708
709
710 /// Parallel distributed vector field which stores the configurational force
711 /// for each fem node corresponding to linear shape function generator (see
712 /// equations 52-53 in
713 /// (https://link.aps.org/doi/10.1103/PhysRevB.97.165132)). This vector
714 /// doesn't contain contribution from terms which have sums over k points.
716
717 /// Parallel distributed vector field which stores the configurational force
718 /// for each fem node corresponding to linear shape function generator (see
719 /// equations 52-53 in
720 /// (https://link.aps.org/doi/10.1103/PhysRevB.97.165132)). This vector only
721 /// containts contribution from the electrostatic part.
723
724#ifdef USE_COMPLEX
725 /// Parallel distributed vector field which stores the configurational force
726 /// for each fem node corresponding to linear shape function generator (see
727 /// equations 52-53 in
728 /// (https://link.aps.org/doi/10.1103/PhysRevB.97.165132)). This vector only
729 /// contains contribution from terms which have sums over k points.
730 distributedCPUVec<double> d_configForceVectorLinFEKPoints;
731#endif
732
733
734 std::vector<double> d_forceAtomsFloating;
735
736#ifdef USE_COMPLEX
737 std::vector<double> d_forceAtomsFloatingKPoints;
738#endif
739
740
741
742 /// Gaussian generator constant. Gaussian generator: Gamma(r)=
743 /// exp(-d_gaussianConstant*r^2)
744 /// FIXME: Until the hanging nodes surface integral issue is fixed use a
745 /// value >=4.0
746 // double d_gaussianConstant;
747
748 /// Storage for configurational force on all global atoms.
749 std::vector<double> d_globalAtomsForces;
750
751 /// Storage for configurational stress tensor
752 dealii::Tensor<2, 3, double> d_stress;
753
754
755
756 /* Part of the stress tensor which is summed over k points.
757 * It is a temporary data structure required for stress evaluation
758 * (d_stress) when parallization over k points is on.
759 */
760 dealii::Tensor<2, 3, double> d_stressKPoints;
761
762 /* Dont use true except for debugging forces only without mesh movement, as
763 * gaussian ovelap on atoms for move mesh is by default set to false
764 */
766
767 /// pointer to dft class
769
770 /// Finite element object for configurational force computation. Linear
771 /// finite elements with three force field components are used.
772 dealii::FESystem<3> FEForce;
773
774 /* DofHandler on which we define the configurational force field. Each
775 * geometric fem node has three dofs corresponding the the three force
776 * components. The generator for the configurational force on the fem node
777 * is the linear shape function attached to it. This DofHandler is based on
778 * the same triangulation on which we solve the dft problem.
779 */
780 dealii::DoFHandler<3> d_dofHandlerForce;
781
782 /* DofHandler on which we define the configurational force field from
783 * electrostatic part (without psp). Each geometric fem node has three dofs
784 * corresponding the the three force components. The generator for the
785 * configurational force on the fem node is the linear shape function
786 * attached to it. This DofHandler is based on the same triangulation on
787 * which we solve the dft problem.
788 */
789 dealii::DoFHandler<3> d_dofHandlerForceElectro;
790
791 /// Index of the d_dofHandlerForce in the MatrixFree object stored in
792 /// dftClass. This is required to correctly use FEEvaluation class.
794
795 /// Index of the d_dofHandlerForceElectro in the MatrixFree object stored in
796 /// dftClass. This is required to correctly use FEEvaluation class.
798
799 /// dealii::IndexSet of locally owned dofs of in d_dofHandlerForce the
800 /// current processor
802
803 /// dealii::IndexSet of locally owned dofs of in d_dofHandlerForceElectro
804 /// the current processor
806
807 /// dealii::IndexSet of locally relevant dofs of in d_dofHandlerForce the
808 /// current processor
810
811 /// dealii::IndexSet of locally relevant dofs of in d_dofHandlerForceElectro
812 /// the current processor
814
815 /// Constraint matrix for hanging node and periodic constaints on
816 /// d_dofHandlerForce.
817 dealii::AffineConstraints<double> d_constraintsNoneForce;
818
819 /// Constraint matrix for hanging node and periodic constaints on
820 /// d_dofHandlerForceElectro.
821 dealii::AffineConstraints<double> d_constraintsNoneForceElectro;
822
823 /// Internal data: map < <atomId,force component>, globaldof in
824 /// d_dofHandlerForce>
825 std::map<std::pair<dftfe::uInt, dftfe::uInt>, dftfe::uInt> d_atomsForceDofs;
826
827 /// Internal data: map < <atomId,force component>, globaldof in
828 /// d_dofHandlerForceElectro>
829 std::map<std::pair<dftfe::uInt, dftfe::uInt>, dftfe::uInt>
831
832 /// Internal data: stores cell iterators of all cells in
833 /// dftPtr->d_dofHandler which are part of the vself ball. Outer vector is
834 /// over vself bins.
835 std::vector<std::vector<dealii::DoFHandler<3>::active_cell_iterator>>
837
838 /// Internal data: stores cell iterators of all cells in d_dofHandlerForce
839 /// which are part of the vself ball. Outer vector is over vself bins.
840 std::vector<std::vector<dealii::DoFHandler<3>::active_cell_iterator>>
842
843 /// Internal data: stores map of vself ball cell Id to the closest atom Id
844 /// of that cell. Outer vector over vself bins.
845 std::vector<std::map<dealii::CellId, dftfe::uInt>>
847
848 /// Internal data: stores the map of atom Id (only in the local processor)
849 /// to the vself bin Id.
850 std::map<dftfe::uInt, dftfe::uInt> d_AtomIdBinIdLocalDofHandler;
851
852 /* Internal data: stores the face ids of dftPtr->d_dofHandler (single
853 * component field) on which to evaluate the vself ball surface integral in
854 * the configurational force expression. Outer vector is over the vself
855 * bins. Inner map is between the cell iterator and the vector of face ids
856 * to integrate on for that cell iterator.
857 */
858 std::vector<std::map<dealii::DoFHandler<3>::active_cell_iterator,
859 std::vector<dftfe::uInt>>>
861
862 /* Internal data: stores the face ids of d_dofHandlerForce (three component
863 * field) on which to evaluate the vself ball surface integral in the
864 * configurational force expression. Outer vector is over the vself bins.
865 * Inner map is between the cell iterator and the vector of face ids to
866 * integrate on for that cell iterator.
867 */
868 std::vector<std::map<dealii::DoFHandler<3>::active_cell_iterator,
869 std::vector<dftfe::uInt>>>
871
872 /// Internal data: stores cell iterators of all cells in
873 /// dftPtr->d_dofHandler which are part of the vself ball. Outer vector is
874 /// over vself bins.
875 std::vector<std::vector<dealii::DoFHandler<3>::active_cell_iterator>>
877
878 /// Internal data: stores cell iterators of all cells in d_dofHandlerForce
879 /// which are part of the vself ball. Outer vector is over vself bins.
880 std::vector<std::vector<dealii::DoFHandler<3>::active_cell_iterator>>
882
883 /// Internal data: stores map of vself ball cell Id to the closest atom Id
884 /// of that cell. Outer vector over vself bins.
885 std::vector<std::map<dealii::CellId, dftfe::uInt>>
887
888 /// Internal data: stores the map of atom Id (only in the local processor)
889 /// to the vself bin Id.
890 std::map<dftfe::uInt, dftfe::uInt> d_AtomIdBinIdLocalDofHandlerElectro;
891
892 /* Internal data: stores the face ids of dftPtr->d_dofHandler (single
893 * component field) on which to evaluate the vself ball surface integral in
894 * the configurational force expression. Outer vector is over the vself
895 * bins. Inner map is between the cell iterator and the vector of face ids
896 * to integrate on for that cell iterator.
897 */
898 std::vector<std::map<dealii::DoFHandler<3>::active_cell_iterator,
899 std::vector<dftfe::uInt>>>
901
902 /* Internal data: stores the face ids of d_dofHandlerForce (three component
903 * field) on which to evaluate the vself ball surface integral in the
904 * configurational force expression. Outer vector is over the vself bins.
905 * Inner map is between the cell iterator and the vector of face ids to
906 * integrate on for that cell iterator.
907 */
908 std::vector<std::map<dealii::DoFHandler<3>::active_cell_iterator,
909 std::vector<dftfe::uInt>>>
911
912 std::map<dealii::CellId, dealii::DoFHandler<3>::active_cell_iterator>
914
915 std::vector<distributedCPUVec<double>> d_gaussianWeightsVecAtoms;
916
917
918
919 /// domain decomposition mpi_communicator
920 const MPI_Comm d_mpiCommParent;
921
922 /// domain decomposition mpi_communicator
923 const MPI_Comm mpi_communicator;
924
926
927 /// number of mpi processes in the current pool
929
930 /// current mpi process id in the current pool
932
933 /// conditional stream object to enable printing only on root processor
934 /// across all pools
935 dealii::ConditionalOStream pcout;
936 };
937
938} // namespace dftfe
939#endif
Definition AtomicCenteredNonLocalOperator.h:58
Definition feevaluationWrapper.h:104
This class is the primary interface location of all other parts of the DFT-FE code for all steps invo...
Definition dft.h:111
Namespace which declares the input parameters and the functions to parse them from the input paramete...
Definition dftParameters.h:36
Calculates dispersion correction to energy, force and stress.
Definition dftd.h:37
void computeConfigurationalForceEEshelbyTensorFPSPFnlLinFE(const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt eigenDofHandlerIndex, const dftfe::uInt smearedChargeQuadratureId, const dftfe::uInt lpspQuadratureIdElectro, const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt phiTotDofHandlerIndexElectro, const distributedCPUVec< double > &phiTotRhoOutElectro, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &rhoOutValues, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &gradRhoOutValues, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValuesLpsp, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValuesLpsp, const std::map< dealii::CellId, std::vector< double > > &rhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &gradRhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &hessianRhoCoreValues, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms, const std::map< dealii::CellId, std::vector< double > > &pseudoVLocElectro, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtomsElectro, const vselfBinsManager &vselfBinsManagerElectro)
void computeConfigurationalForceEEshelbyEElectroPhiTot(const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt phiTotDofHandlerIndexElectro, const dftfe::uInt smearedChargeQuadratureId, const dftfe::uInt lpspQuadratureIdElectro, const distributedCPUVec< double > &phiTotRhoOutElectro, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValues, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValuesLpsp, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValues, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValuesLpsp, const std::map< dealii::CellId, std::vector< double > > &pseudoVLocElectro, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtomsElectro, const vselfBinsManager &vselfBinsManagerElectro)
std::vector< std::map< dealii::CellId, dftfe::uInt > > d_cellsVselfBallsClosestAtomIdDofHandler
Definition force.h:846
void computeConfigurationalForceTotalLinFE(const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt eigenDofHandlerIndex, const dftfe::uInt smearedChargeQuadratureId, const dftfe::uInt lpspQuadratureIdElectro, const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt phiTotDofHandlerIndexElectro, const distributedCPUVec< double > &phiTotRhoOutElectro, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &rhoOutValues, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &gradRhoOutValues, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValuesLpsp, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValuesLpsp, const std::map< dealii::CellId, std::vector< double > > &rhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &gradRhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &hessianRhoCoreValues, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms, const std::map< dealii::CellId, std::vector< double > > &pseudoVLocElectro, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtomsElectro, const vselfBinsManager &vselfBinsManagerElectro)
void computeAtomsForcesGaussianGenerator(bool allowGaussianOverlapOnAtoms=false)
void stressEnlElementalContribution(dealii::Tensor< 2, 3, double > &stressContribution, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt numQuadPoints, const std::vector< double > &jxwQuadsSubCells, const dftfe::uInt cell, const dftfe::uInt numNonLocalAtomsCurrentProcess, const std::shared_ptr< AtomicCenteredNonLocalOperator< dataTypes::number, memorySpace > > nonLocalOp, const std::vector< dftfe::uInt > &numberPseudoWaveFunctionsPerAtom, const std::map< dealii::CellId, dftfe::uInt > &cellIdToCellNumberMap, const std::vector< dataTypes::number > &zetalmDeltaVlProductDistImageAtoms, const std::vector< dataTypes::number > &projectorKetTimesPsiTimesVTimesPartOccContractionGradPsiQuadsFlattened, const bool isSpinPolarized)
std::vector< std::map< dealii::DoFHandler< 3 >::active_cell_iterator, std::vector< dftfe::uInt > > > d_cellFacesVselfBallSurfacesDofHandlerForceElectro
Definition force.h:910
std::vector< std::vector< dealii::DoFHandler< 3 >::active_cell_iterator > > d_cellsVselfBallsDofHandlerForce
Definition force.h:841
void configForceLinFEInit(const dealii::MatrixFree< 3, double > &matrixFreeData, const dealii::MatrixFree< 3, double > &matrixFreeDataElectro)
void computeStressEEshelbyEPSPEnlEk(const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt eigenDofHandlerIndex, const dftfe::uInt smearedChargeQuadratureId, const dftfe::uInt lpspQuadratureIdElectro, const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt phiTotDofHandlerIndexElectro, const distributedCPUVec< double > &phiTotRhoOutElectro, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &rhoOutValues, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &gradRhoOutValues, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValuesLpsp, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValuesLpsp, const std::map< dealii::CellId, std::vector< double > > &pseudoVLocElectro, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtomsElectro, const std::map< dealii::CellId, std::vector< double > > &rhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &gradRhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &hessianRhoCoreValues, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms, const vselfBinsManager &vselfBinsManagerElectro)
std::vector< std::map< dealii::DoFHandler< 3 >::active_cell_iterator, std::vector< dftfe::uInt > > > d_cellFacesVselfBallSurfacesDofHandler
Definition force.h:860
void addEVselfSmearedStressContribution(FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &gradVselfQuads, const std::vector< dftfe::uInt > &nonTrivialAtomImageIdsMacroCell, const std::map< dealii::CellId, std::vector< dftfe::Int > > &bQuadAtomIdsAllAtomsImages, const dealii::AlignedVector< dealii::VectorizedArray< double > > &smearedbQuads)
void distributeForceContributionFPSPLocalGammaAtoms(const std::map< dftfe::uInt, std::vector< double > > &forceContributionFPSPLocalGammaAtoms, const std::map< std::pair< dftfe::uInt, dftfe::uInt >, dftfe::uInt > &atomsForceDofs, const dealii::AffineConstraints< double > &constraintsNoneForce, distributedCPUVec< double > &configForceVectorLinFE)
distributedCPUVec< double > d_configForceVectorLinFE
Definition force.h:715
void createBinObjectsForce(const dealii::DoFHandler< 3 > &dofHandler, const dealii::DoFHandler< 3 > &dofHandlerForce, const dealii::AffineConstraints< double > &hangingPlusPBCConstraints, const vselfBinsManager &vselfBinsManager, std::vector< std::vector< dealii::DoFHandler< 3 >::active_cell_iterator > > &cellsVselfBallsDofHandler, std::vector< std::vector< dealii::DoFHandler< 3 >::active_cell_iterator > > &cellsVselfBallsDofHandlerForce, std::vector< std::map< dealii::CellId, dftfe::uInt > > &cellsVselfBallsClosestAtomIdDofHandler, std::map< dftfe::uInt, dftfe::uInt > &AtomIdBinIdLocalDofHandler, std::vector< std::map< dealii::DoFHandler< 3 >::active_cell_iterator, std::vector< dftfe::uInt > > > &cellFacesVselfBallSurfacesDofHandler, std::vector< std::map< dealii::DoFHandler< 3 >::active_cell_iterator, std::vector< dftfe::uInt > > > &cellFacesVselfBallSurfacesDofHandlerForce)
void initMoved(std::vector< const dealii::DoFHandler< 3 > * > &dofHandlerVectorMatrixFree, std::vector< const dealii::AffineConstraints< double > * > &constraintsVectorMatrixFree, const bool isElectrostaticsMesh)
initializes data structures inside forceClass which depend on the moved mesh.
const MPI_Comm d_mpiCommParent
domain decomposition mpi_communicator
Definition force.h:920
void computeConfigurationalForceEselfNoSurfaceLinFE()
dealii::IndexSet d_locally_owned_dofsForce
Definition force.h:801
std::vector< std::map< dealii::CellId, dftfe::uInt > > d_cellsVselfBallsClosestAtomIdDofHandlerElectro
Definition force.h:886
void printStress()
prints the currently stored configurational stress tensor.
void computeConfigurationalForceEselfLinFE(const dealii::DoFHandler< 3 > &dofHandlerElectro, const vselfBinsManager &vselfBinsManagerElectro, const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt smearedChargeQuadratureId)
dealii::IndexSet d_locally_relevant_dofsForceElectro
Definition force.h:813
std::map< dftfe::uInt, dftfe::uInt > d_AtomIdBinIdLocalDofHandlerElectro
Definition force.h:890
dealii::DoFHandler< 3 > d_dofHandlerForce
Definition force.h:780
void initUnmoved(const dealii::Triangulation< 3, 3 > &triangulation, const dealii::Triangulation< 3, 3 > &serialTriangulation, const std::vector< std::vector< double > > &domainBoundingVectors, const bool isElectrostaticsMesh)
initializes data structures inside forceClass assuming unmoved triangulation.
dftClass< memorySpace > * dftPtr
pointer to dft class
Definition force.h:768
void addENonlinearCoreCorrectionStressContribution(FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::VectorizedArray< double > > &vxcQuads, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &derExcGradRho, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms)
void FPhiTotSmearedChargesGammaAtomsElementalContribution(std::map< dftfe::uInt, std::vector< double > > &forceContributionSmearedChargesGammaAtoms, FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &gradPhiTotQuads, const std::vector< dftfe::uInt > &nonTrivialAtomIdsMacroCell, const std::map< dealii::CellId, std::vector< dftfe::Int > > &bQuadAtomIdsAllAtoms, const dealii::AlignedVector< dealii::VectorizedArray< double > > &smearedbQuads)
const dftfe::uInt n_mpi_processes
number of mpi processes in the current pool
Definition force.h:928
std::vector< double > d_globalAtomsForces
Storage for configurational force on all global atoms.
Definition force.h:749
std::vector< std::vector< dealii::DoFHandler< 3 >::active_cell_iterator > > d_cellsVselfBallsDofHandler
Definition force.h:836
dealii::Tensor< 2, 3, double > & getStress()
returns a copy of the current stress tensor value.
void accumulateForceContributionGammaAtomsFloating(const std::map< dftfe::uInt, std::vector< double > > &forceContributionLocalGammaAtoms, std::vector< double > &accumForcesVector)
std::vector< double > d_forceAtomsFloating
Definition force.h:734
std::map< std::pair< dftfe::uInt, dftfe::uInt >, dftfe::uInt > d_atomsForceDofsElectro
Definition force.h:830
void addEPhiTotSmearedStressContribution(FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &gradPhiTotQuads, const std::vector< dftfe::uInt > &nonTrivialAtomImageIdsMacroCell, const std::map< dealii::CellId, std::vector< dftfe::Int > > &bQuadAtomIdsAllAtomsImages, const dealii::AlignedVector< dealii::VectorizedArray< double > > &smearedbQuads)
const dftParameters & d_dftParams
Definition force.h:925
void initPseudoData()
initializes and precomputes pseudopotential related data structuers required for configurational forc...
std::vector< std::map< dealii::DoFHandler< 3 >::active_cell_iterator, std::vector< dftfe::uInt > > > d_cellFacesVselfBallSurfacesDofHandlerElectro
Definition force.h:900
dealii::FESystem< 3 > FEForce
Definition force.h:772
void computeStressEEshelbyEElectroPhiTot(const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt phiTotDofHandlerIndexElectro, const dftfe::uInt smearedChargeQuadratureId, const dftfe::uInt lpspQuadratureIdElectro, const distributedCPUVec< double > &phiTotRhoOutElectro, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValues, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValuesLpsp, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValuesElectro, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValuesElectroLpsp, const std::map< dealii::CellId, std::vector< double > > &pseudoVLocElectro, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtomsElectro, const vselfBinsManager &vselfBinsManagerElectro)
void FNonlinearCoreCorrectionGammaAtomsElementalContribution(std::map< dftfe::uInt, std::vector< double > > &forceContributionFNonlinearCoreCorrectionGammaAtoms, FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &derExcGradRho, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms)
void locateAtomCoreNodesForce(const dealii::DoFHandler< 3 > &dofHandlerForce, const dealii::IndexSet &locally_owned_dofsForce, std::map< std::pair< dftfe::uInt, dftfe::uInt >, dftfe::uInt > &atomsForceDofs)
get the value of Gaussian generator parameter (d_gaussianConstant). Gaussian generator: Gamma(r)= exp...
void distributeForceContributionFnlGammaAtoms(const std::map< dftfe::uInt, std::vector< double > > &forceContributionFnlGammaAtoms)
forceClass(dftClass< memorySpace > *_dftPtr, const MPI_Comm &mpi_comm_parent, const MPI_Comm &mpi_comm_domain, const dftParameters &dftParams)
Constructor.
void computeConfigurationalForcePhiExtLinFE()
dealii::Tensor< 2, 3, double > d_stress
Storage for configurational stress tensor.
Definition force.h:752
void configForceLinFEFinalize()
void computeFloatingAtomsForces()
dftfe::uInt d_forceDofHandlerIndexElectro
Definition force.h:797
dftfe::uInt d_forceDofHandlerIndex
Definition force.h:793
void addEPSPStressContribution(dealii::FEValues< 3 > &feValues, dealii::FEFaceValues< 3 > &feFaceValues, FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt phiTotDofHandlerIndexElectro, const dftfe::uInt cell, const dealii::AlignedVector< dealii::VectorizedArray< double > > &rhoQuads, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &gradRhoQuads, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtoms, const vselfBinsManager &vselfBinsManager, const std::vector< std::map< dealii::CellId, dftfe::uInt > > &cellsVselfBallsClosestAtomIdDofHandler)
std::vector< distributedCPUVec< double > > d_gaussianWeightsVecAtoms
Definition force.h:915
void computeElementalNonLocalPseudoOVDataForce()
distributedCPUVec< double > d_configForceVectorLinFEElectro
Definition force.h:722
void FNonlinearCoreCorrectionGammaAtomsElementalContribution(std::map< dftfe::uInt, std::vector< double > > &forceContributionFNonlinearCoreCorrectionGammaAtoms, FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::VectorizedArray< double > > &vxcQuads, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms)
std::map< std::pair< dftfe::uInt, dftfe::uInt >, dftfe::uInt > d_atomsForceDofs
Definition force.h:825
const dftfe::uInt this_mpi_process
current mpi process id in the current pool
Definition force.h:931
std::map< dealii::CellId, dealii::DoFHandler< 3 >::active_cell_iterator > d_cellIdToActiveCellIteratorMapDofHandlerRhoNodalElectro
Definition force.h:913
const bool d_allowGaussianOverlapOnAtoms
Definition force.h:765
void FnlGammaAtomsElementalContribution(std::map< dftfe::uInt, std::vector< double > > &forceContributionFnlGammaAtoms, const dealii::MatrixFree< 3, double > &matrixFreeData, FEEvaluationWrapperClass< 3 > &forceEvalNLP, const std::shared_ptr< AtomicCenteredNonLocalOperator< dataTypes::number, memorySpace > > nonLocalOp, dftfe::uInt numNonLocalAtomsCurrentProcess, const std::vector< dftfe::Int > &globalChargeIdNonLocalAtoms, const std::vector< dftfe::uInt > &numberPseudoWaveFunctionsPerAtom, const dftfe::uInt cell, const std::map< dealii::CellId, dftfe::uInt > &cellIdToCellNumberMap, const std::vector< dataTypes::number > &zetaDeltaVQuadsFlattened, const std::vector< dataTypes::number > &projectorKetTimesPsiTimesVTimesPartOccContractionGradPsiQuadsFlattened)
dealii::ConditionalOStream pcout
Definition force.h:935
void FnlGammaxElementalContribution(dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &FVectQuads, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt numQuadPoints, const std::shared_ptr< AtomicCenteredNonLocalOperator< dataTypes::number, memorySpace > > nonLocalOp, const dftfe::uInt numNonLocalAtomsCurrentProcess, const std::vector< dftfe::Int > &globalChargeIdNonLocalAtoms, const std::vector< dftfe::uInt > &numberPseudoWaveFunctionsPerAtom, const dftfe::uInt cell, const std::map< dealii::CellId, dftfe::uInt > &cellIdToCellNumberMap, const std::vector< dataTypes::number > &zetaDeltaVQuadsFlattened, const std::vector< dataTypes::number > &projectorKetTimesPsiTimesVTimesPartOccContractionGradPsiQuadsFlattened)
void FPSPLocalGammaAtomsElementalContribution(std::map< dftfe::uInt, std::vector< double > > &forceContributionFPSPLocalGammaAtoms, dealii::FEValues< 3 > &feValues, dealii::FEFaceValues< 3 > &feFaceValues, FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt phiTotDofHandlerIndexElectro, const dftfe::uInt cell, const dealii::AlignedVector< dealii::VectorizedArray< double > > &rhoQuads, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &gradRhoQuads, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtoms, const vselfBinsManager &vselfBinsManager, const std::vector< std::map< dealii::CellId, dftfe::uInt > > &cellsVselfBallsClosestAtomIdDofHandler)
dealii::AffineConstraints< double > d_constraintsNoneForce
Definition force.h:817
dealii::IndexSet d_locally_owned_dofsForceElectro
Definition force.h:805
dealii::IndexSet d_locally_relevant_dofsForce
Definition force.h:809
void printAtomsForces()
prints the currently stored configurational forces on atoms and the Gaussian generator constant used ...
void computeAtomsForces(const dealii::MatrixFree< 3, double > &matrixFreeData, const dispersionCorrection &dispersionCorr, const dftfe::uInt eigenDofHandlerIndex, const dftfe::uInt smearedChargeQuadratureId, const dftfe::uInt lpspQuadratureIdElectro, const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt phiTotDofHandlerIndexElectro, const distributedCPUVec< double > &phiTotRhoOutElectro, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &rhoOutValues, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &gradRhoOutValues, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValuesLpsp, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValuesLpsp, const std::map< dealii::CellId, std::vector< double > > &rhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &gradRhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &hessianRhoCoreValues, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms, const std::map< dealii::CellId, std::vector< double > > &pseudoVLocElectro, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtomsElectro, const dealii::AffineConstraints< double > &hangingPlusPBCConstraintsElectro, const vselfBinsManager &vselfBinsManagerElectro)
computes the configurational force on all atoms corresponding to a Gaussian generator,...
std::vector< std::map< dealii::DoFHandler< 3 >::active_cell_iterator, std::vector< dftfe::uInt > > > d_cellFacesVselfBallSurfacesDofHandlerForce
Definition force.h:870
std::map< dftfe::uInt, dftfe::uInt > d_AtomIdBinIdLocalDofHandler
Definition force.h:850
std::vector< std::vector< dealii::DoFHandler< 3 >::active_cell_iterator > > d_cellsVselfBallsDofHandlerElectro
Definition force.h:876
std::vector< std::vector< dealii::DoFHandler< 3 >::active_cell_iterator > > d_cellsVselfBallsDofHandlerForceElectro
Definition force.h:881
const MPI_Comm mpi_communicator
domain decomposition mpi_communicator
Definition force.h:923
std::vector< double > & getAtomsForces()
returns a copy of the configurational force on all global atoms.
dealii::Tensor< 2, 3, double > d_stressKPoints
Definition force.h:760
void addENonlinearCoreCorrectionStressContributionSpinPolarized(FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::VectorizedArray< double > > &vxcQuadsSpin0, const dealii::AlignedVector< dealii::VectorizedArray< double > > &vxcQuadsSpin1, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &derExcGradRhoSpin0, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &derExcGradRhoSpin1, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms, const bool isXCGGA=false)
dealii::AffineConstraints< double > d_constraintsNoneForceElectro
Definition force.h:821
void computeStressEself(const dealii::DoFHandler< 3 > &dofHandlerElectro, const vselfBinsManager &vselfBinsManagerElectro, const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt smearedChargeQuadratureId)
void FNonlinearCoreCorrectionGammaAtomsElementalContributionSpinPolarized(std::map< dftfe::uInt, std::vector< double > > &forceContributionFNonlinearCoreCorrectionGammaAtoms, FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::VectorizedArray< double > > &vxcQuadsSpin0, const dealii::AlignedVector< dealii::VectorizedArray< double > > &vxcQuadsSpin1, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &derExcGradRhoSpin0, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &derExcGradRhoSpin1, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms, const bool isXCGGA=false)
void computeStress(const dealii::MatrixFree< 3, double > &matrixFreeData, const dispersionCorrection &dispersionCorr, const dftfe::uInt eigenDofHandlerIndex, const dftfe::uInt smearedChargeQuadratureId, const dftfe::uInt lpspQuadratureIdElectro, const dealii::MatrixFree< 3, double > &matrixFreeDataElectro, const dftfe::uInt phiTotDofHandlerIndexElectro, const distributedCPUVec< double > &phiTotRhoOutElectro, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &rhoOutValues, const std::vector< dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > > &gradRhoOutValues, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &rhoTotalOutValuesLpsp, const dftfe::utils::MemoryStorage< double, dftfe::utils::MemorySpace::HOST > &gradRhoTotalOutValuesLpsp, const std::map< dealii::CellId, std::vector< double > > &pseudoVLocElectro, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &pseudoVLocAtomsElectro, const std::map< dealii::CellId, std::vector< double > > &rhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &gradRhoCoreValues, const std::map< dealii::CellId, std::vector< double > > &hessianRhoCoreValues, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &gradRhoCoreAtoms, const std::map< dftfe::uInt, std::map< dealii::CellId, std::vector< double > > > &hessianRhoCoreAtoms, const dealii::AffineConstraints< double > &hangingPlusPBCConstraintsElectro, const vselfBinsManager &vselfBinsManagerElectro)
Update force generator Gaussian constant.
void FVselfSmearedChargesGammaAtomsElementalContribution(std::map< dftfe::uInt, std::vector< double > > &forceContributionSmearedChargesGammaAtoms, FEEvaluationWrapperClass< 3 > &forceEval, const dealii::MatrixFree< 3, double > &matrixFreeData, const dftfe::uInt cell, const dealii::AlignedVector< dealii::Tensor< 1, 3, dealii::VectorizedArray< double > > > &gradVselfBinQuads, const std::vector< dftfe::uInt > &nonTrivialAtomIdsMacroCell, const std::map< dealii::CellId, std::vector< dftfe::Int > > &bQuadAtomIdsAllAtoms, const dealii::AlignedVector< dealii::VectorizedArray< double > > &smearedbQuads)
dealii::DoFHandler< 3 > d_dofHandlerForceElectro
Definition force.h:789
Definition MemoryStorage.h:33
Categorizes atoms into bins for efficient solution of nuclear electrostatic self-potential.
Definition vselfBinsManager.h:34
Definition pseudoPotentialToDftfeConverter.cc:34
dealii::LinearAlgebra::distributed::Vector< elem_type, dealii::MemorySpace::Host > distributedCPUVec
Definition headers.h:92
std::uint32_t uInt
Definition TypeConfig.h:10