Bridge++  Ver. 1.1.x
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
test_Spectrum_CRSMatrix_Clover_Lexical.cpp
Go to the documentation of this file.
1 
14 #include "parameterManager_YAML.h"
15 
16 #include "bridgeIO.h"
17 using Bridge::vout;
18 
19 #include "fopr_CRS.h"
20 #include "fopr_Clover.h"
21 
22 #include "source_4spinor_Local.h"
23 
24 #include "solver_CG.h"
25 
26 #include "gaugeConfig.h"
27 
28 #ifdef USE_TEST
29 #include "test.h"
30 #endif
31 
32 #ifdef USE_TESTMANAGER_AUTOREGISTER
33 #include "testManager.h"
34 #endif
35 
36 //====================================================================
38 
53 namespace Test_Spectrum_CRSMatrix {
54  //- test-private parameters
55  namespace {
56  const std::string filename_input = "test_Spectrum_CRSMatrix_Clover_Lexical.yaml";
57  const std::string filename_output = "stdout";
58 
59  class Parameters_Test_Spectrum_CRSMatrix : public Parameters {
60  public:
61  Parameters_Test_Spectrum_CRSMatrix()
62  {
63  Register_string("gauge_config_type_input", "NULL");
64  Register_string("config_filename_input", "NULL");
65 
66  Register_string("matrix_output", "NULL");
67  Register_string("source_output", "NULL");
68  Register_string("solution_output", "NULL");
69 
70  Register_string("verbose_level", "NULL");
71 
72  Register_double("expected_result", 0.0);
73  }
74  };
75  }
76 
77  //- prototype declaration
78  int clover_lex(void);
79 
80  int CRSsolver(
81  const string& solution,
82  const string& matrix,
83  const string& source,
84  double& result /* return value */
85  );
86 
87 
88 #ifdef USE_TESTMANAGER_AUTOREGISTER
89 #ifdef USE_MPI
90  // this test runs only in single-node environment.
91 #else
92  //- NB. CRS test is skipped, because it is time-consuming.
93  // namespace {
94  // static const bool is_registered = TestManager::RegisterTest(
95  // "Spectrum.CRSMatrix.Clover_Lexical",
96  // clover_lex
97  // );
98  // }
99 #endif
100 #endif
101 
102  //====================================================================
103  int clover_lex(void)
104  {
105  // #### parameter setup ####
106  int Nc = CommonParameters::Nc();
107  int Nd = CommonParameters::Nd();
108  int Ndim = CommonParameters::Ndim();
109  int Nvol = CommonParameters::Nvol();
110 
111  Parameters_Test_Spectrum_CRSMatrix params_test;
112 
113  Parameters params_all;
114 
115  params_all.Register_Parameters("Test_Spectrum_CRSMatrix", &params_test);
116 
117  ParameterManager_YAML params_manager;
118  params_manager.read_params(filename_input, &params_all);
119 
120  const string str_gconf_read = params_test.get_string("gauge_config_type_input");
121  const string readfile = params_test.get_string("config_filename_input");
122  const string matrix_file = params_test.get_string("matrix_output");
123  const string source_file = params_test.get_string("source_output");
124  const string solution_file = params_test.get_string("solution_output");
125  const string str_vlevel = params_test.get_string("verbose_level");
126 #ifdef USE_TEST
127  const double expected_result = params_test.get_double("expected_result");
128 #endif
129 
130  Bridge::VerboseLevel vl = vout.set_verbose_level(str_vlevel);
131 
132  //- print input parameters
133  vout.general(vl, " gconf_read = %s\n", str_gconf_read.c_str());
134  vout.general(vl, " readfile = %s\n", readfile.c_str());
135  vout.general(vl, " matrix_output = %s\n", matrix_file.c_str());
136  vout.general(vl, " source_output = %s\n", source_file.c_str());
137  vout.general(vl, " solution_output = %s\n", solution_file.c_str());
138  vout.general(vl, " vlevel = %s\n", str_vlevel.c_str());
139  vout.general(vl, "\n");
140 
141  //- input parameter check
142  int err = 0;
143  err += ParameterCheck::non_NULL(str_gconf_read);
144  err += ParameterCheck::non_NULL(readfile);
145  err += ParameterCheck::non_NULL(matrix_file);
146  err += ParameterCheck::non_NULL(source_file);
147  err += ParameterCheck::non_NULL(solution_file);
148 
149  if (err) {
150  vout.crucial(vl, "Test_Spectrum_CRSMatrix: Input parameters have not been set.\n");
151  abort();
152  }
153 
154 
155  // #### Set up a gauge configuration ####
156  Field_G *U = new Field_G(Nvol, Ndim);
157  GaugeConfig *gconf_read = new GaugeConfig(str_gconf_read);
158  gconf_read->read_file((Field *)U, readfile);
159  // gconf_read->set_cold((Field*)U);
160 
161 
162  // #### object setup #####
163  Fopr_Clover *fopr_c = new Fopr_Clover();
164  double kappa = 0.12;
165  double cSW = 1.0;
166  valarray<int> boundary(Ndim);
167  boundary[0] = -1;
168  boundary[1] = -1;
169  boundary[2] = -1;
170  boundary[3] = -1;
171  fopr_c->set_parameters(kappa, cSW, boundary);
172  fopr_c->set_config(U);
173  fopr_c->set_mode("D");
174 
175  Fopr_CRS *fopr_crs = new Fopr_CRS(fopr_c);
176  fopr_crs->write_matrix(matrix_file);
177 
178  int Niter = 5000;
179  double Stop_cond = 1.0e-28;
180  Solver *solver = new Solver_CG((Fopr *)fopr_c);
181  solver->set_parameters(Niter, Stop_cond);
182 
183  valarray<int> source_position(Ndim);
184  source_position[0] = 0;
185  source_position[1] = 0;
186  source_position[2] = 0;
187  source_position[3] = 0;
188 
189  Source_4spinor_Local *source = new Source_4spinor_Local();
190  source->set_parameters(source_position);
191 
192 
193  // #### Execution main part ####
194  valarray<Field_F> sq(Nc * Nd);
195  Field_F xq, b, b2;
196 
197  {
198  int ispin = 0;
199  {
200  int icolor = 0;
201  //for(int ispin = 0; ispin < Nd; ++ispin){
202  // for(int icolor = 0; icolor < Nc; ++icolor){
203 
204  source->set(b, icolor, ispin);
205 
206  b.write_text("source4x8_clover.crs");
207 
208  fopr_c->set_mode("D");
209  b2 = (Field_F)fopr_c->mult_dag(b);
210 
211  int Nconv;
212  double diff_CG;
213  fopr_c->set_mode("DdagD");
214  solver->solve(xq, b2, Nconv, diff_CG);
215  vout.general(vl, " ispin = %2d icolor = %2d Nconv = %4d diff = %12.6e\n",
216  ispin, icolor, Nconv, diff_CG);
217 
218  xq.write_text("solution4x8_clover.crs");
219 
220  Field y(b);
221  fopr_c->set_mode("D");
222  y -= fopr_c->mult(xq);
223  double yy = y.norm2();
224  vout.general(vl, " standard norm2 = %.8e\n", yy);
225 
226  sq[icolor + Nc * ispin] = xq;
227  }
228  }
229 
230 #if 0
231  // ### Solver for CRS version ###
232 
233  Solver *solver_crs = new Solver_CG((Fopr *)fopr_crs2);
234  solver_crs->set_parameters(Niter, Stop_cond);
235 
236  double yy = 0.0L;
237  {
238  int ispin = 0;
239  {
240  int icolor = 0;
241  //for(int ispin = 0; ispin < Nd; ++ispin){
242  // for(int icolor = 0; icolor < Nc; ++icolor){
243 
244  source->set(b, icolor, ispin);
245 
246  fopr_crs2->set_mode("Ddag");
247  b2 = (Field_F)fopr_crs2->mult((Field)b);
248 
249  int Nconv;
250  double diff_CG;
251  fopr_crs2->set_mode("DdagD");
252  solver_crs->solve(xq, b2, Nconv, diff_CG);
253  vout.general(vl, " ispin = %2d icolor = %2d Nconv = %4d diff = %12.6e\n",
254  ispin, icolor, Nconv, diff_CG);
255 
256  Field_F y(b);
257  fopr_crs2->set_mode("D");
258  y -= (Field_F)fopr_crs2->mult((Field)xq);
259  yy = y.norm2();
260  vout.general(vl, " standard norm2 = %.8e\n", yy);
261 
262  sq[icolor + Nc * ispin] = xq;
263  }
264  }
265 #endif
266 
267 
268  double result;
269  CRSsolver(solution_file, matrix_file, source_file, result);
270 
271 
272  // #### tidy up ####
273  delete U;
274  delete gconf_read;
275 
276  delete source;
277  delete solver;
278 
279  delete fopr_crs;
280  delete fopr_c;
281 
282 #if 0
283  delete solver_crs;
284  delete fopr_crs2;
285 #endif
286 
287 
288 #ifdef USE_TEST
289  return Test::verify(expected_result, result);
290 
291 #else
292  return EXIT_SUCCESS;
293 #endif
294  }
295 } // namespace Test_Spectrum_CRSMatrix