00001
00008 #include <math.h>
00009 #include "SourcesAndLosses.h"
00010
00011 #define err 1e-3
00012
00013 SourcesAndLosses::SourcesAndLosses(ParamStructure::SL parameters, Grid &grid, int numberOfIterations, double timeStep) {
00014 Initialize(parameters, grid, numberOfIterations, timeStep);
00015 }
00016
00024 void SourcesAndLosses::Initialize(ParamStructure::SL parameters, Grid &grid, int numberOfIterations, double timeStep) {
00025 SL_parameters = parameters;
00026
00027 int it;
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044 int il, im, ia;
00045
00046 for (it = 0; it < numberOfIterations; it++) {
00047
00048 (*this).push_back(Matrix3D<double>());
00049
00050 (*this)[it].AllocateMemory(grid.L.size, grid.epc.size, grid.alpha.size);
00051 for (il = 0; il < grid.L.size; il++) {
00052 for (im = 0; im < grid.epc.size; im++) {
00053 for (ia = 0; ia < grid.alpha.size; ia++) {
00054 (*this)[it][il][im][ia] = 0;
00055 }
00056 }
00057 }
00058 }
00059
00060
00061 if(parameters.SL_L_top) {
00062
00063 Matrix3D<double> temp(numberOfIterations, grid.pc.size, grid.alpha.size);
00064 Matrix3D<double> loaded_time(numberOfIterations, grid.pc.size, grid.alpha.size);
00065 Matrix3D<double> loaded_epc(numberOfIterations, grid.pc.size, grid.alpha.size);
00066 Matrix3D<double> loaded_alpha(numberOfIterations, grid.pc.size, grid.alpha.size);
00067
00069 temp.readFromFile(parameters.SL_L_top_filename, loaded_time, loaded_epc, loaded_alpha);
00070
00071
00072 for (it = 1; it < numberOfIterations; it++) {
00073 for (im = 0; im < grid.epc.size; im++) {
00074 for (ia = 0; ia < grid.alpha.size; ia++) {
00075
00076 if (fabs(timeStep*it - loaded_time[it][im][ia]) > 1
00077 || fabs(grid.epc[grid.L.size-1][im][ia] - loaded_epc[it][im][ia]) > err
00078 || fabs(grid.alpha[grid.L.size-1][im][ia] - loaded_alpha[it][im][ia]) > err) {
00079 throw error_msg("L_TOP_LOAD_GRID_MISMATCH", "Loading %s: grid mismatch [%d, %d, %d].\nLoaded: %e, %e, %e\nGrid: %e, %e, %e\n",
00080 parameters.SL_L_top_filename.c_str(), it, im, ia,
00081 loaded_time[it][im][ia], loaded_epc[it][im][ia], loaded_alpha[it][im][ia],
00082 (timeStep*it), grid.epc[grid.L.size-1][im][ia], grid.alpha[grid.L.size-1][im][ia]);
00083 } else {
00084
00085
00086
00087 (*this)[it][grid.L.size-2][im][ia] = (*this)[it][grid.L.size-2][im][ia] + temp[it][im][ia];
00088 }
00089 }
00090 }
00091 }
00092 }
00093
00094
00095 if(parameters.SL_E_min) {
00096
00097 Matrix3D<double> temp(numberOfIterations, grid.L.size, grid.alpha.size);
00098 Matrix3D<double> loaded_time(numberOfIterations, grid.L.size, grid.alpha.size);
00099 Matrix3D<double> loaded_L(numberOfIterations, grid.L.size, grid.alpha.size);
00100 Matrix3D<double> loaded_alpha(numberOfIterations, grid.L.size, grid.alpha.size);
00101
00103 temp.readFromFile(parameters.SL_E_min_filename, loaded_time, loaded_L, loaded_alpha);
00104
00105 for (it = 1; it < numberOfIterations; it++) {
00106 for (il = 0; il < grid.L.size; il++) {
00107 for (ia = 0; ia < grid.alpha.size; ia++) {
00108
00109 if (fabs(timeStep*it - loaded_time[it][il][ia]) > 1
00110 || fabs(grid.L[il][1][ia] - loaded_L[it][il][ia]) > err
00111 || fabs(grid.alpha[il][1][ia] - loaded_alpha[it][il][ia]) > err) {
00112 throw error_msg("L_TOP_LOAD_GRID_MISMATCH", "Loading %s: grid mismatch [%d, %d, %d].\nLoaded: %e, %e, %e\nGrid: %e, %e, %e\n",
00113 parameters.SL_E_min_filename.c_str(), it, il, ia,
00114 loaded_time[it][il][ia], loaded_L[it][il][ia], loaded_alpha[it][il][ia],
00115 (timeStep*it), grid.L[il][1][ia], grid.alpha[il][1][ia]);
00116 } else {
00117 (*this)[it][il][1][ia] = (*this)[it][il][1][ia] + temp[it][il][ia];
00118 }
00119 }
00120 }
00121 }
00122 }
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132 }