VERB_code_2.3
BoundaryConditions.cpp
Go to the documentation of this file.
1 
8 #include "BoundaryConditions.h"
9 #include <math.h>
10 
11 // Function 'strncasecmp' called '_strnicmp' in windows
12 #if defined(_WIN32) || defined(_WIN64)
13 #define strncasecmp _strnicmp
14 #define strcasecmp _stricmp
15 #endif
16 
17 
18 
30  int ix, iy;
31 
32  if (!arr.initialized) {
33  // Allocate memory
34  arr.AllocateMemory(psd2DSlice.size_x, psd2DSlice.size_y);
35  }
36 
37  // store type in a class variable, for the ease of future use
38  this->type = parameters.type;
39  if (parameters.type == "BCT_CONSTANT_VALUE") {
40  // Saving boundary value inside the class.
41  // If boundary is on constant value of the function - saving this type and the value.
42  this->calculationType = "BCT_CONSTANT_VALUE";
43 
44  // If value is zero, we don't put zero to PSD, but put some small value
45  //(*this) = VF::max(value, VC::zero_f);
46  arr = VF::max(parameters.value, VC::zero_f);
47 
48  } else if (parameters.type == "BCT_CONSTANT_PSD") {
49  // If boundary is on constant PSD value - saving this type and we need to store the initial values of PSD later, after initialization of PSD.
50  this->calculationType = "BCT_CONSTANT_VALUE";
51 
52  // So, this line store psd slice as a boundry, and change all zeros to some small value
53  for (ix = 0; ix < arr.size_x; ix++) {
54  for (iy = 0; iy < arr.size_y; iy++) {
55  arr[ix][iy] = VF::max(psd2DSlice[ix][iy], VC::zero_f);
56  }
57  }
58  } else if (parameters.type == "BCT_CONSTANT_DERIVATIVE") {
59  // If boundary is on constant derivative of the function - saving this type and the value of the derivative (usually - zero derivative).
60  this->calculationType = "BCT_CONSTANT_DERIVATIVE";
61 
62  // boundary condition is a constant derivative - so, storing derivative value for boundary conditions for each boundary point
63  // storing it inside the parent class Matrix.
64  // this - is a pointer to current class. '=' - operator.
65  arr = parameters.value;
66  } else if (parameters.type == "BCT_CONSTANT_SECOND_DERIVATIVE") {
67  // If boundary is on constant second derivative of the function
68  this->calculationType = "BCT_CONSTANT_SECOND_DERIVATIVE";
69 
70  // boundary condition is a constant derivative - so, storing derivative value for boundary conditions for each boundary point
71  // storing it inside the parent class Matrix.
72  // this - is a pointer to current class. '=' - operator.
73  arr = parameters.value;
74  } else if (parameters.type == "BCT_FILE" || parameters.type == "BCT_FILE_GRID" || parameters.type == "BCT_MULTIPLE_FILES") {
75  // If boundary are loaded from file... see BoundaryCondition::LoadBoundaryCondition
76  this->calculationType = "BCT_CONSTANT_VALUE";
77  // this->filename = parameters.filename; // don't need this i guess... but it does not do any bad.
78  } else {
79  throw error_msg("BCT_ERROR", "Unknown boundary condition type: %d", parameters.type.c_str());
80  }
81 }
82 
83 
94  int ix, iy;
95  double err = 1.e-3;
96  string inBuf;
97 
98  if (!arr.initialized) {
99  // Allocating memory
100  arr.AllocateMemory(gridElement1.size_x, gridElement1.size_y);
101  }
102 
103  Output::echo(0, "Loading boundary conditions from files... ");
104  // store type in a class variable, for the ease of future use
105  this->type = parameters.type;
106  //if (initialType == "BCT_FILE" || initialType == "BCT_FILE_GRID") {
107  if (parameters.type == "BCT_FILE" || parameters.type == "BCT_FILE_GRID") {
108  this->calculationType = "BCT_CONSTANT_VALUE";
109  // this->readFromFile(filename);
110  // making stream
111  {
112  ifstream input(parameters.filename.c_str());
113  // checking if file was opened
114  if (input == NULL) {
115  throw error_msg("BC_FILE_OPEN_ERR", "Boundary conditions file %s not found.", parameters.filename.c_str());
116  }
117  // Skip header!
118  input >> inBuf;
119  while (strcasecmp(inBuf.c_str(), "ZONE") != 0 && !input.eof() ) input >> inBuf;
120  input.ignore(9999, '\n'); // read to the end of the line with 'zone'
121  double loadedTime, loadedAx1, loadedAx2, loadedValue;
122  // size_x - time steps, size_y and size_z - two dimensions of PSD slice
123  for (ix = 0; ix < arr.size_x; ix++) {
124  for (iy = 0; iy < arr.size_y; iy++) {
125  if (parameters.type == "BCT_FILE_GRID") input >> loadedAx1 >> loadedAx2;
126  input >> loadedValue;
127  // check if grid is the same
128  if (parameters.type == "BCT_FILE_GRID" && (fabs(loadedAx1 - gridElement1[ix][iy]) > err || fabs(loadedAx2 - gridElement2[ix][iy]) > err)) {
129  if (input.eof()) {
130  throw error_msg("BC_LOAD_GRID_ERR", "Loading %s: not enough data points in the boundary condition file.\n", parameters.filename.c_str());
131  } else {
132  throw error_msg("BC_LOAD_GRID_ERR", "Loading %s: grid mismatch.\nLoaded: ax1=%e, ax2=%e\nGrid: ax1=%e, ax2=%e\n(i1=%d, i2=%d)\n", parameters.filename.c_str(), loadedAx1, loadedAx2, gridElement1[ix][iy], gridElement2[ix][iy], ix, iy);
133  }
134  } else {
135  arr[ix][iy] = loadedValue;
136  }
137  // skip till the end of the line
138  input.ignore(9999, '\n');
139  }
140  }
141  input.close();
142  }
143  arr.writeToFile("./Debug/" + arr.name + "_BC.dat");
144  Output::echo(0, "done.\n");
145  } else if (parameters.type == "BCT_MULTIPLE_FILES") {
146  // In that case we need to change boundary conditions (BC) with time.
147  // Now we load the list of boundary conditions (files with boundary values) and corresponding times, when we need to apply the BC
148  this->calculationType = "BCT_CONSTANT_VALUE";
149 
150  Output::echo(0, "Loading a list of boundary conditions... ");
151  ifstream input(parameters.filename.c_str());
152  // checking if file was opened
153  if (input == NULL) {
154  throw error_msg("BC_FILE_OPEN_ERR", "Boundary conditions file %s not found.", parameters.filename.c_str());
155  }
156 
157  // Skip header of the file
158  input >> inBuf;
159  while (strcasecmp(inBuf.c_str(), "ZONE") != 0 && !input.eof() ) input >> inBuf;
160  input.ignore(9999, '\n'); // read to the end of the line with 'zone'
161 
162  // reading the file
163  double time_tmp;
164  string filename_tmp;
165  while (!input.eof()) {
166  // reading time values and function values from the file into vectors
167  input >> time_tmp >> filename_tmp;
168  BC_change_time.push_back(time_tmp);
169  BC_filename.push_back(filename_tmp);
170  }
171 
172  // Load the first on
173  ifstream input2(BC_filename[0].c_str());
174  // checking if file was opened
175  if (input2 == NULL) {
176  throw error_msg("BC_FILE_OPEN_ERR", "Boundary conditions file %s not found.", parameters.filename.c_str());
177  }
178  Output::echo(3, "Updating a boundary condition: %s", BC_filename[0].c_str());
179  // Skip header of the file
180  string inBuf;
181  input2 >> inBuf;
182  while (strcasecmp(inBuf.c_str(), "ZONE") != 0 && !input2.eof() ) input2 >> inBuf;
183  input2.ignore(9999, '\n'); // read to the end of the line with 'zone'
184  double tmp;
185  for (ix = 0; ix < arr.size_x; ix++) {
186  for (iy = 0; iy < arr.size_y; iy++) {
187  input2 >> arr[ix][iy];
188  }
189  }
190 
191  }
192 
193 }
202 void BoundaryCondition::Update(int iteration, Matrix2D<double> PSD_2D_Slice, double time, double dt) {
203  int ix, iy;
204 
205  if (type == "BCT_MULTIPLE_FILES") {
206  if (time == -1 || dt == -1)
207  throw error_msg("BC_ERR", "To use BCT_MULTIPLE_FILES, time and dt need to be specified when the function is called.");
208  // Find the next time to change BC
209  int i;
210  for (i = 0; i < BC_change_time.size() && BC_change_time[i] < time; i++) {};
211  // Check if it's time to load next BC, i.e. the time to change is between the previous and current time steps
212  if (time - dt < BC_change_time[i] && time >= BC_change_time[i]) {
213  // Load next BC
214  ifstream input(BC_filename[i].c_str());
215  // checking if file was opened
216  if (input == NULL) {
217  throw error_msg("BC_FILE_OPEN_ERR", "Boundary conditions file %s not found.", BC_filename[i].c_str());
218  }
219  Output::echo(3, "Updating a boundary condition: %s", BC_filename[i].c_str());
220  // Skip header of the file
221  string inBuf;
222  input >> inBuf;
223  while (strcasecmp(inBuf.c_str(), "ZONE") != 0 && !input.eof() ) input >> inBuf;
224  input.ignore(9999, '\n'); // read to the end of the line with 'zone'
225  double tmp;
226  for (ix = 0; ix < arr.size_x; ix++) {
227  for (iy = 0; iy < arr.size_y; iy++) {
228  input >> tmp;
229  arr[ix][iy] = tmp;
230  }
231  }
232  }
233  }
234 
235 }
bool initialized
Flag, equal true if initialized.
Definition: Matrix.h:138
int size_y
size x, size y
Definition: Matrix.h:141
double max(double v1, double v2)
Return maximum.
vector< string > BC_filename
If we change boundary conditions with time, this is the filename for the new boundary condition value...
void AllocateMemory(int size_x, int size_y)
Definition: Matrix.cpp:524
string calculationType
Type: constant function/constant derivative.
string filename
File name to load boundary conditions from.
Definition: Parameters.h:177
Matrix2D< double > arr
2D matrix of boundary conditions
int size_x
size x, size y
Definition: Matrix.h:141
void echo(int msglvl, const char *format,...)
Basic function! Write the message to the file.
Definition: Output.cpp:44
Boundary conditions parameters structure.
Definition: Parameters.h:171
string type
Initial type: how to set up boundary conditions. (Like read from a file, etc)
vector< double > BC_change_time
If we change boundary conditions with time, this is the time when we need to change it...
void MakeBoundaryCondition(Parameters_structure::BoundaryCondition parameters, Matrix2D< double > psd2DSlice, Matrix2D< double > gridElement1, Matrix2D< double > gridElement2)
string type
Type. Check StrToVal(string input, BoundaryConditionTypes &place) for known values.
Definition: Parameters.h:173
double value
Value (if it is constant value on the boundary).
Definition: Parameters.h:175
Parameters_structure parameters
Parameters structure, with all parameters from the parameters.ini file. The default parameters define...
Definition: Main.cpp:187
Error message - stack of single_errors.
Definition: error.h:54
string name
name of the Matrix
Definition: Matrix.h:143
void writeToFile(string filename)
Definition: Matrix.cpp:682
static const double zero_f
Minimum value for PSD.
void Update(int iteration, Matrix2D< double > PSD_2D_Slice, double time=-1, double dt=-1)
void LoadBoundaryCondition(Parameters_structure::BoundaryCondition parameters, Matrix2D< double > gridElement1, Matrix2D< double > gridElement2)