PSCF v1.2
rpc/field/FieldIo.tpp
1#ifndef RPC_FIELD_IO_TPP
2#define RPC_FIELD_IO_TPP
3
4/*
5* PSCF - Polymer Self-Consistent Field Theory
6*
7* Copyright 2016 - 2022, The Regents of the University of Minnesota
8* Distributed under the terms of the GNU General Public License.
9*/
10
11#include "FieldIo.h"
12#include <prdc/field/FieldIoReal.tpp> // base class implementation
13
14#include <rpc/field/Domain.h>
15
16#include <prdc/field/fieldIoUtil.h>
17#include <prdc/crystal/fieldHeader.h>
18#include <prdc/crystal/Basis.h>
19#include <prdc/crystal/UnitCell.h>
20#include <prdc/cpu/complex.h>
21#include <prdc/cpu/types.h>
22#include <pscf/mesh/Mesh.h>
23#include <pscf/math/IntVec.h>
24
25namespace Pscf {
26namespace Rpc {
27
28 using namespace Util;
29 using namespace Pscf::Prdc;
30 using namespace Pscf::Prdc::Cpu;
31
32 /*
33 * Constructor.
34 */
35 template <int D>
37 : FieldIoReal< D, RField<D>, RFieldDft<D>, FFT<D> >()
38 {}
39
40 /*
41 * Destructor.
42 */
43 template <int D>
46
47 /*
48 * Read an array of fields in r-grid format.
49 */
50 template <int D>
52 std::istream &in,
53 DArray<RField<D> >& fields,
54 UnitCell<D>& unitCell) const
55 {
56 // Read header
57 int nMonomer;
58 bool isSymmetric;
59 readFieldHeader(in, nMonomer, unitCell, isSymmetric);
60 readMeshDimensions(in, mesh().dimensions());
61 checkAllocateFields(fields, nMonomer, mesh().dimensions());
62
63 // Read data
64 // Rpg:: Allocate host arrays
65 Prdc::readRGridData(in, fields, nMonomer, mesh().dimensions());
66 // Rpg:: Copy host -> device
67
68 }
69
70 /*
71 * Read the data section of an array of fields in r-grid format.
72 */
73 template <int D>
75 std::istream& in,
76 DArray< RField<D> >& fields,
77 int nMonomer) const
78 {
79 checkAllocateFields(fields, nMonomer, mesh().dimensions());
80 // Rpg:: Allocate host arrays
81 Prdc::readRGridData(in, fields, nMonomer, mesh().dimensions());
82 // Rpg:: Copy host -> device
83 }
84
85 /*
86 * Read a single field in r-grid format.
87 */
88 template <int D>
90 std::istream &in,
91 RField<D> & field,
92 UnitCell<D>& unitCell) const
93 {
94
95 // Read header
96 int nMonomer;
97 bool isSymmetric;
98 readFieldHeader(in, nMonomer, unitCell, isSymmetric);
99 UTIL_CHECK(nMonomer == 1);
100 readMeshDimensions(in, mesh().dimensions());
101
102 // Read data
103 // Rpg:: Allocate host arrays
104 checkAllocateField(field, mesh().dimensions());
105 Prdc::readRGridData(in, field, mesh().dimensions());
106 // Rpg:: Copy host -> device
107 }
108
109 /*
110 * Write an array of fields in r-grid format.
111 */
112 template <int D>
114 std::ostream &out,
115 DArray<RField<D> > const & fields,
116 UnitCell<D> const & unitCell,
117 bool writeHeader,
118 bool isSymmetric,
119 bool writeMeshSize) const
120 {
121 // Inspect fields array, infer nMonomer and meshDimensions
122 int nMonomer;
123 IntVec<D> meshDimensions;
124 inspectFields(fields, nMonomer, meshDimensions);
125
126 // Write header
127 if (writeHeader){
128 writeFieldHeader(out, nMonomer, unitCell, isSymmetric);
129 }
130 if (writeMeshSize){
131 writeMeshDimensions(out, meshDimensions);
132 }
133
134 // Write data section
135 // Rpg:: Allocate host arrays
136 // Rpg:: Copy device -> host
137 Prdc::writeRGridData(out, fields, nMonomer, meshDimensions);
138 }
139
140 /*
141 * Write a single field in r-grid format.
142 */
143 template <int D>
145 std::ostream &out,
146 RField<D> const & field,
147 UnitCell<D> const & unitCell,
148 bool writeHeader,
149 bool isSymmetric) const
150 {
151 IntVec<D> meshDimensions = field.meshDimensions();
152
153 // Write header
154 if (writeHeader) {
155 writeFieldHeader(out, 1, unitCell, isSymmetric);
156 writeMeshDimensions(out, meshDimensions);
157 }
158
159 // Write data
160 // Rpg:: Allocate host array
161 // Rpg:: Copy device -> host
162 Prdc::writeRGridData(out, field, meshDimensions);
163 }
164
165 /*
166 * Read an array of fields in k-grid format
167 */
168 template <int D>
170 std::istream &in,
171 DArray<RFieldDft<D> >& fields,
172 UnitCell<D>& unitCell) const
173 {
174 // Read header
175 int nMonomer;
176 bool isSymmetric;
177 readFieldHeader(in, nMonomer, unitCell, isSymmetric);
178 readMeshDimensions(in, mesh().dimensions());
179
180 checkAllocateFields(fields, nMonomer, mesh().dimensions());
181 IntVec<D> dftDimensions = fields[0].dftDimensions();
182
183 // Read data
184 // Rpg:: Allocate host arrays
185 Prdc::readKGridData(in, fields, nMonomer, dftDimensions);
186 // Rpg:: Copy host -> device
187 }
188
189 /*
190 * Write an array of fields in k-grid format
191 */
192 template <int D>
194 std::ostream &out,
195 DArray<RFieldDft<D> > const & fields,
196 UnitCell<D> const & unitCell,
197 bool isSymmetric) const
198 {
199 // Inspect fields array - infer nMonomer and dimensions
200 int nMonomer;
201 IntVec<D> meshDimensions;
202 inspectFields(fields, nMonomer, meshDimensions);
203 IntVec<D> dftDimensions = fields[0].dftDimensions();
204
205 // Write file
206 writeFieldHeader(out, nMonomer, unitCell, isSymmetric);
207 writeMeshDimensions(out, meshDimensions);
208
209 // Write data
210 // Rpg:: Allocate host arrays
211 // Rpg:: Copy device -> host
212 Prdc::writeKGridData(out, fields, nMonomer, dftDimensions);
213 }
214
215 /*
216 * Convert an array of fields from basis to k-grid format.
217 */
218 template <int D>
220 DArray<double> const & in,
221 RFieldDft<D>& out) const
222 {
223 // Rpg: Allocate host array
224 Prdc::convertBasisToKGrid(in, out, basis(), out.dftDimensions());
225 // Rpg: Copy host -> device
226 }
227
228 /*
229 * Convert an array of fields from k-grid to basis format.
230 */
231 template <int D>
233 RFieldDft<D> const & in,
234 DArray<double>& out,
235 bool checkSymmetry,
236 double epsilon) const
237 {
238 // Rpg: Allocate host arrays
239 // Rpg: Copy device -> host
240 Prdc::convertKGridToBasis(in, out, basis(), in.dftDimensions(),
241 checkSymmetry, epsilon);
242 }
243
244 /*
245 * Test if an real field DFT has the declared space group symmetry.
246 */
247 template <int D>
249 RFieldDft<D> const & in,
250 double epsilon,
251 bool verbose) const
252 {
253 // Rpg:: Allocate host array
254 // Rpg: Copy device -> host
255 return Prdc::hasSymmetry(in, basis(), in.dftDimensions(),
256 epsilon, verbose);
257 }
258
259 /*
260 * Test if an real field DFT has the declared space group symmetry.
261 */
262 template <int D>
264 DArray<double> & field,
265 double factor) const
266 {
267 int n = field.capacity();
268 for (int i = 0; i < n; ++i) {
269 field[i] *= factor;
270 }
271 }
272
273 /*
274 * Test if an real field DFT has the declared space group symmetry.
275 */
276 template <int D>
278 RField<D> & field,
279 double factor) const
280 {
281 int n = field.capacity();
282 for (int i = 0; i < n; ++i) {
283 field[i] *= factor;
284 }
285 }
286
287 /*
288 * Replicate the unit cell for an array of r-grid fields.
289 */
290 template <int D>
292 std::ostream &out,
293 DArray< RField<D> > const & fields,
294 UnitCell<D> const & unitCell,
295 IntVec<D> const & replicas) const
296
297 {
298 // Inspect fields to obtain nMonomer and meshDimensions
299 int nMonomer;
300 IntVec<D> meshDimensions;
301 inspectFields(fields, nMonomer, meshDimensions);
302
303 // Rpg: Allocate hostArrays
304 // Rpg: Copy device -> host
305 Prdc::replicateUnitCell(out, fields, meshDimensions,
306 unitCell, replicas);
307 }
308
309 /*
310 * Expand spatial dimension of an array of r-grid fields.
311 */
312 template <int D>
314 std::ostream &out,
315 DArray< RField<D> > const & fields,
316 UnitCell<D> const & unitCell,
317 int d,
318 DArray<int> const& newGridDimensions) const
319 {
320 IntVec<D> meshDimensions = fields[0].meshDimensions();
321
322 // Rpg: Allocate hostArrays
323 // Rpg: Copy device -> host
324 Prdc::expandRGridDimension(out, fields, meshDimensions,
325 unitCell, d, newGridDimensions);
326 }
327
328}
329}
330#endif
An IntVec<D, T> is a D-component vector of elements of integer type T.
Definition IntVec.h:27
Fourier transform wrapper.
Fourier transform of a real field on an FFT mesh.
IntVec< D > const & dftDimensions() const
Return vector of dft (Fourier) grid dimensions by constant reference.
Field of real double precision values on an FFT mesh.
const IntVec< D > & meshDimensions() const
Return mesh dimensions by constant reference.
Definition cpu/RField.h:112
File input/output operations and format conversions for fields.
Definition FieldIoReal.h:68
Base template for UnitCell<D> classes, D=1, 2 or 3.
Definition rpg/System.h:34
void expandRGridDimension(std::ostream &out, DArray< RField< D > > const &fields, UnitCell< D > const &unitCell, int d, DArray< int > const &newGridDimensions) const override
Expand spatial dimension of an array of r-grid fields.
void replicateUnitCell(std::ostream &out, DArray< RField< D > > const &fields, UnitCell< D > const &unitCell, IntVec< D > const &replicas) const override
Write r-grid fields in a replicated unit cell to std::ostream.
bool hasSymmetry(RFieldDft< D > const &in, double epsilon=1.0e-8, bool verbose=true) const override
Check if a k-grid field has the declared space group symmetry.
void scaleFieldRGrid(RField< D > &field, double factor) const override
Rescale a single r-grid field by a scalar factor.
void writeFieldsKGrid(std::ostream &out, DArray< RFieldDft< D > > const &fields, UnitCell< D > const &unitCell, bool isSymmetric=true) const override
Write array of RFieldDft objects (k-space fields) to file.
void readFieldsRGrid(std::istream &in, DArray< RField< D > > &fields, UnitCell< D > &unitCell) const override
Read array of RField objects (r-grid fields) from a stream.
void convertBasisToKGrid(DArray< double > const &components, RFieldDft< D > &dft) const override
Convert a field from symmetrized basis to Fourier grid (k-grid).
void readFieldsRGridData(std::istream &in, DArray< RField< D > > &fields, int nMonomer) const override
Read data for an array of r-grid fields, with no header section.
void readFieldsKGrid(std::istream &in, DArray< RFieldDft< D > > &fields, UnitCell< D > &unitCell) const override
Read array of RFieldDft objects (k-space fields) from a stream.
void convertKGridToBasis(RFieldDft< D > const &in, DArray< double > &out, bool checkSymmetry=true, double epsilon=1.0e-8) const override
Convert a field from Fourier (k-grid) to symmetrized basis form.
void writeFieldsRGrid(std::ostream &out, DArray< RField< D > > const &fields, UnitCell< D > const &unitCell, bool writeHeader=true, bool isSymmetric=true, bool writeMeshSize=true) const override
Write array of RField objects (fields on r-space grid) to a stream.
void scaleFieldBasis(DArray< double > &field, double factor) const override
Rescale a single field in basis format by a scalar factor.
void readFieldRGrid(std::istream &in, RField< D > &field, UnitCell< D > &unitCell) const override
Read a single RField (field on an r-space grid) from a stream.
FieldIo()
Default constructor.
void writeFieldRGrid(std::ostream &out, RField< D > const &field, UnitCell< D > const &unitCell, bool writeHeader=true, bool isSymmetric=true) const override
Write a single RField (field on an r-space grid) to a stream.
int capacity() const
Return allocated size.
Definition Array.h:159
Dynamically allocatable contiguous array template.
#define UTIL_CHECK(condition)
Assertion macro suitable for serial or parallel production code.
Definition global.h:68
void writeFieldHeader(std::ostream &out, int ver1, int ver2, UnitCell< D > const &cell, std::string const &groupName, int nMonomer)
Write common part of field header (fortran PSCF format).
void readFieldHeader(std::istream &in, int &ver1, int &ver2, UnitCell< D > &cell, std::string &groupName, int &nMonomer)
Read common part of field header (fortran PSCF format).
void replicateUnitCell(std::ostream &out, DArray< AT > const &fields, IntVec< D > const &meshDimensions, UnitCell< D > const &unitCell, IntVec< D > const &replicas)
Write r-grid fields in a replicated unit cell to std::ostream.
void writeKGridData(std::ostream &out, DArray< AT > const &fields, int nMonomer, IntVec< D > const &dftDimensions)
Write data for array of k-grid fields, with no header section.
bool hasSymmetry(AT const &in, Basis< D > const &basis, IntVec< D > const &dftDimensions, double epsilon=1.0e-8, bool verbose=true)
Check if a k-grid field has the declared space group symmetry.
void inspectFields(DArray< FT > const &fields, int &nMonomer, IntVec< D > &dimensions)
Inspect dimensions of a DArray of fields, each of type FT.
void readKGridData(std::istream &in, DArray< AT > &fields, int nMonomer, IntVec< D > const &dftDimensions)
Read data for array of k-grid fields, with no header section.
void convertBasisToKGrid(DArray< double > const &components, AT &dft, Basis< D > const &basis, IntVec< D > const &dftDimensions)
Convert a real field from symmetrized basis to Fourier grid.
void writeMeshDimensions(std::ostream &out, IntVec< D > const &meshDimensions)
Write mesh dimensions to a field file header.
void expandRGridDimension(std::ostream &out, DArray< AT > const &fields, IntVec< D > const &meshDimensions, UnitCell< D > const &unitCell, int d, DArray< int > newGridDimensions)
Expand the dimensionality of space from D to d.
void writeRGridData(std::ostream &out, DArray< AT > const &fields, int nMonomer, IntVec< D > const &dimensions)
Write data for array of r-grid fields, with no header section.
void convertKGridToBasis(AT const &in, DArray< double > &out, Basis< D > const &basis, IntVec< D > const &dftDimensions, bool checkSymmetry=true, double epsilon=1.0e-8)
Convert a real field from Fourier grid to symmetrized basis.
void readRGridData(std::istream &in, DArray< AT > &fields, int nMonomer, IntVec< D > const &dimensions)
Read data for array of r-grid fields, with no header section.
void readMeshDimensions(std::istream &in, IntVec< D > const &meshDimensions)
Read mesh dimensions from a field file header.
void checkAllocateFields(DArray< FT > &fields, int nMonomer, IntVec< D > const &dimensions)
Check allocation of an array of fields, allocate if necessary.
void checkAllocateField(FT &field, IntVec< D > const &dimensions)
Check allocation of a single field, allocate if necessary.
Fields and FFTs for periodic boundary conditions (CPU)
Definition CField.cpp:12
Periodic fields and crystallography.
Definition CField.cpp:11
PSCF package top-level namespace.
Definition param_pc.dox:1
Utility classes for scientific computation.