1#ifndef RPG_LR_POST_AM_COMPRESSOR_TPP
2#define RPG_LR_POST_AM_COMPRESSOR_TPP
11#include "LrAmCompressor.h"
12#include <rpg/System.h>
13#include <rpg/fts/compressor/intra/IntraCorrelation.h>
14#include <prdc/cuda/resources.h>
15#include <pscf/mesh/MeshIterator.h>
28 isIntraCalculated_(false),
30 { setClassName(
"LrAmCompressor"); }
50 const int nMonomer = system().mixture().nMonomer();
51 const int meshSize = system().domain().mesh().size();
52 IntVec<D> const & dimensions = system().mesh().dimensions();
58 for (
int i = 0; i < D; ++i) {
60 kMeshDimensions_[i] = dimensions[i];
62 kMeshDimensions_[i] = dimensions[i]/2 + 1;
68 for (
int i = 0; i < D; ++i) {
69 kSize_ *= kMeshDimensions_[i];
74 w0_.allocate(nMonomer);
75 wFieldTmp_.allocate(nMonomer);
76 resid_.allocate(dimensions);
77 residK_.allocate(dimensions);
78 intraCorrelationK_.allocate(kMeshDimensions_);
79 for (
int i = 0; i < nMonomer; ++i) {
80 w0_[i].allocate(dimensions);
81 wFieldTmp_[i].allocate(dimensions);
88 if (!isIntraCalculated_){
89 intra_.computeIntraCorrelations(intraCorrelationK_);
90 isIntraCalculated_ =
true;
94 for (
int i = 0; i < nMonomer; ++i) {
128 double LrAmCompressor<D>::maxAbs(DeviceArray<cudaReal>
const & a)
136 LrAmCompressor<D>::updateBasis(
RingBuffer< DeviceArray<cudaReal> > & basis,
137 RingBuffer< DeviceArray<cudaReal> >
const & hists)
143 if (basis[0].isAllocated()) {
144 UTIL_CHECK(basis[0].capacity() == hists[0].capacity());
146 basis[0].allocate(hists[0].capacity());
154 LrAmCompressor<D>::addHistories(DeviceArray<cudaReal>& trial,
155 RingBuffer<DeviceArray<cudaReal> >
const & basis,
159 for (
int i = 0; i < nHist; i++) {
165 double LrAmCompressor<D>::computeLambda(
double r)
171 void LrAmCompressor<D>::addPredictedError(DeviceArray<cudaReal>& fieldTrial,
172 DeviceArray<cudaReal>
const & resTrial,
175 int n = fieldTrial.capacity();
176 const double vMonomer = system().mixture().vMonomer();
182 system().fft().forwardTransform(resid_, residK_);
188 system().fft().inverseTransformUnsafe(residK_, resid_);
196 bool LrAmCompressor<D>::hasInitialGuess()
197 {
return system().w().hasData(); }
201 int LrAmCompressor<D>::nElements()
202 {
return system().domain().mesh().size(); }
206 void LrAmCompressor<D>::getCurrent(DeviceArray<cudaReal>& curr)
219 void LrAmCompressor<D>::evaluate()
227 void LrAmCompressor<D>::getResidual(DeviceArray<cudaReal>& resid)
229 const int nMonomer = system().mixture().nMonomer();
235 for (
int i = 1; i < nMonomer; i++) {
242 void LrAmCompressor<D>::update(DeviceArray<cudaReal>& newGuess)
245 const int nMonomer = system().mixture().nMonomer();
248 for (
int i = 0; i < nMonomer; i++) {
253 system().setWRGrid(wFieldTmp_);
257 void LrAmCompressor<D>::outputToLog()
265 out <<
"LrAmCompressor time contributions:\n";
Template for Anderson mixing iterator algorithm.
Dynamic array on the GPU device with aligned data.
int capacity() const
Return allocated capacity.
An IntVec<D, T> is a D-component vector of elements of integer type T.
Base class for iterators that impose incompressibility.
Anderson Mixing compressor with linear-response mixing step.
~LrAmCompressor()
Destructor.
void readParameters(std::istream &in)
Read all parameters and initialize.
int compress()
Compress to obtain partial saddle point w+.
void clearTimers()
Clear all timers (reset accumulated time to zero).
void setup(bool isContinuation)
Initialize just before entry to iterative loop.
LrAmCompressor(System< D > &system)
Constructor.
void outputTimers(std::ostream &out)
Return compressor times contributions.
Main class for calculations that represent one system.
Dynamically allocatable contiguous array template.
Class for storing history of previous values in an array.
File containing preprocessor macros for error handling.
#define UTIL_CHECK(condition)
Assertion macro suitable for serial or parallel production code.
cudaReal innerProduct(DeviceArray< cudaReal > const &a, DeviceArray< cudaReal > const &b)
Compute inner product of two real arrays (GPU kernel wrapper).
cudaReal maxAbs(DeviceArray< cudaReal > const &in)
Get maximum absolute magnitude of array elements (GPU kernel wrapper).
void addEqV(DeviceArray< cudaReal > &a, DeviceArray< cudaReal > const &b, const int beginIdA, const int beginIdB, const int n)
Vector addition in-place, a[i] += b[i], kernel wrapper (cudaReal).
void eqV(DeviceArray< cudaReal > &a, DeviceArray< cudaReal > const &b, const int beginIdA, const int beginIdB, const int n)
Vector assignment, a[i] = b[i], kernel wrapper (cudaReal).
void subVS(DeviceArray< cudaReal > &a, DeviceArray< cudaReal > const &b, cudaReal const c, const int beginIdA, const int beginIdB, const int n)
Vector subtraction, a[i] = b[i] - c, kernel wrapper (cudaReal).
void addVV(DeviceArray< cudaReal > &a, DeviceArray< cudaReal > const &b, DeviceArray< cudaReal > const &c, const int beginIdA, const int beginIdB, const int beginIdC, const int n)
Vector addition, a[i] = b[i] + c[i], kernel wrapper (cudaReal).
void subVV(DeviceArray< cudaReal > &a, DeviceArray< cudaReal > const &b, DeviceArray< cudaReal > const &c, const int beginIdA, const int beginIdB, const int beginIdC, const int n)
Vector subtraction, a[i] = b[i] - c[i], kernel wrapper (cudaReal).
void divEqVc(DeviceArray< cudaComplex > &a, DeviceArray< cudaReal > const &b, cudaReal const c)
Vector division in-place w/ coeff., a[i] /= (b[i] * c), kernel wrapper.
void addEqVc(DeviceArray< cudaReal > &a, DeviceArray< cudaReal > const &b, cudaReal const c)
Vector addition in-place w/ coefficient, a[i] += b[i] * c, kernel wrapper.
PSCF package top-level namespace.
Utility classes for scientific computation.