123 #define _CRT_SECURE_NO_DEPRECATE
124 #define __VERB_VERSION_NUMBER__ "2.3.0.0" // Daniel Choi version
144 #include "../Logging/Output.h"
147 #include "../VariousFunctions/variousFunctions.h"
150 #include "../Diffusion/PSD.h"
153 #include "../Grid/Grid.h"
156 #include "../Grid/BoundaryConditions.h"
159 #include "../Grid/AdditionalSourcesAndLosses.h"
162 #include "../Parameters/Parameters.h"
166 #if defined(_WIN32) || defined(_WIN64)
167 #define _CRTDBG_MAP_ALLOC
195 int main(
int argc,
char* argv[]) {
212 radialDiffusionGrid, localDiffusionsGrid;
216 psdRadialDiffusion, psdLocalDiffusions;
224 L_LowerBoundaryCondition, L_UpperBoundaryCondition,
225 pc_LowerBoundaryCondition, pc_UpperBoundaryCondition,
226 alpha_LowerBoundaryCondition, alpha_UpperBoundaryCondition;
241 double simulation_time = 0;
249 std::clock_t zero_time, start_time, current_time;
251 time_t now = time(NULL);
253 printf(
"%d\n", localtime(&now)->tm_year);
255 if (localtime(&now)->tm_year + 1900 >= 2016) {
256 cout <<
"The program version is outdated, please update." << endl;
267 Output::echo(0,
"Compiled time: %s %s\n", __DATE__, __TIME__);
304 parameters.radialDiffusionGrid_pc,
305 parameters.radialDiffusionGrid_alpha,
306 parameters.radialDiffusionGrid_epc,
314 localDiffusionsGrid.
Create_Grid(parameters.localDiffusionsGrid_L,
315 parameters.localDiffusionsGrid_pc,
316 parameters.localDiffusionsGrid_alpha,
317 parameters.localDiffusionsGrid_epc,
367 if (parameters.L_UpperBoundaryCondition.
type ==
"BCT_FILE" || parameters.L_UpperBoundaryCondition.
type ==
"BCT_FILE_GRID" || parameters.L_UpperBoundaryCondition.
type ==
"BCT_MULTIPLE_FILES") {
369 parameters.L_UpperBoundaryCondition,
373 if (parameters.pc_LowerBoundaryCondition.
type ==
"BCT_FILE" || parameters.pc_LowerBoundaryCondition.
type ==
"BCT_FILE_GRID" || parameters.pc_LowerBoundaryCondition.
type ==
"BCT_MULTIPLE_FILES") {
375 parameters.pc_LowerBoundaryCondition,
414 parameters.L_LowerBoundaryCondition,
419 parameters.L_UpperBoundaryCondition,
426 parameters.pc_LowerBoundaryCondition,
431 parameters.pc_UpperBoundaryCondition,
437 parameters.alpha_LowerBoundaryCondition,
442 parameters.alpha_UpperBoundaryCondition,
465 psdLocalDiffusions.
Interpolate(psdRadialDiffusion, parameters.
interpolation, radialDiffusionGrid, localDiffusionsGrid, pc_LowerBoundaryCondition.
arr, pc_UpperBoundaryCondition.
arr);
502 Output1DHeaders(output1D, Daa, Dpcpc, Dpca, DaaLpp, DpcpcLpp, DpcaLpp);
504 Output1DValues(output1D, Daa, Dpcpc, Dpca, DaaLpp, DpcpcLpp, DpcaLpp, simulation_time, parameters, 0);
521 simulation_time = parameters.
timeStep * iteration;
524 Output::echo(1,
"Step %d of %d. Time = %d days, %.0f hours.\n", iteration, parameters.
totalIterationsNumber,
int(simulation_time), (simulation_time-
int(simulation_time))*24);
529 if (radialDiffusionGrid.
L.
size > 1) {
541 L_LowerBoundaryCondition.
Update(iteration, psdRadialDiffusion.
arr.
xSlice(0));
542 L_UpperBoundaryCondition.
Update(iteration, psdRadialDiffusion.
arr.
xSlice(radialDiffusionGrid.
L.
size-1), simulation_time, parameters.
timeStep);
556 Output::echo(3,
"Updating radial diffusion coefficient...\n");
557 DLL.
MakeDLL(radialDiffusionGrid.
L, radialDiffusionGrid.
pc, radialDiffusionGrid.
alpha, parameters.
Kp[iteration], parameters.
DLLType);
583 radialDiffusionGrid.
L, radialDiffusionGrid.
pc, radialDiffusionGrid.
alpha, radialDiffusionGrid.
Jacobian,
584 L_LowerBoundaryCondition.
arr,
585 L_UpperBoundaryCondition.
arr * parameters.
Bf[iteration],
588 parameters.
tau[iteration],
589 parameters.
tauLpp[iteration]);
620 psdLocalDiffusions.
Interpolate(psdRadialDiffusion, parameters.
interpolation, radialDiffusionGrid, localDiffusionsGrid, pc_LowerBoundaryCondition.
arr, psdRadialDiffusion.
arr.
ySlice(radialDiffusionGrid.
pc.
size-1));
656 if (localDiffusionsGrid.
pc.
size > 1) {
658 pc_LowerBoundaryCondition.
Update(iteration, psdLocalDiffusions.
arr.
ySlice(0));
659 pc_UpperBoundaryCondition.
Update(iteration, psdLocalDiffusions.
arr.
ySlice(radialDiffusionGrid.
pc.
size-1));
662 if (localDiffusionsGrid.
alpha.
size > 1) {
664 alpha_LowerBoundaryCondition.
Update(iteration, psdLocalDiffusions.
arr.
zSlice(0));
672 Output::echo(3,
"Updating momentum diffusion coefficients...\n");
678 Output::echo(3,
"Updating pitch-angle diffusion coefficients...\n");
684 Output::echo(3,
"Updating mixed diffusion coefficients... \n");
713 if (localDiffusionsGrid.
pc.
size > 1 && localDiffusionsGrid.
alpha.
size > 1) {
714 Output::echo(3,
"Calculating local diffusions using BLOCK method... \n");
718 parameters.
Lpp[iteration],
722 localDiffusionsGrid.
L, localDiffusionsGrid.
pc, localDiffusionsGrid.
alpha, localDiffusionsGrid.
Jacobian,
723 pc_LowerBoundaryCondition.
arr,
724 pc_UpperBoundaryCondition.
arr,
725 alpha_LowerBoundaryCondition.
arr,
726 alpha_UpperBoundaryCondition.
arr,
768 Output::echo(3,
"Calculating energy diffusions using SPLIT method... \n");
770 parameters.
Lpp[iteration],
772 localDiffusionsGrid.
L, localDiffusionsGrid.
pc, localDiffusionsGrid.
alpha, localDiffusionsGrid.
Jacobian,
773 pc_LowerBoundaryCondition.
arr,
774 pc_UpperBoundaryCondition.
arr,
785 Output::echo(3,
"Calculating pitch angle diffusions using SPLIT method... \n");
787 parameters.
Lpp[iteration],
789 localDiffusionsGrid.
L, localDiffusionsGrid.
pc, localDiffusionsGrid.
alpha, localDiffusionsGrid.
Jacobian,
790 alpha_LowerBoundaryCondition.
arr,
791 alpha_UpperBoundaryCondition.
arr,
802 Output::echo(3,
"Calculating mixed terms EXPLICITLY ... \n");
804 parameters.
Lpp[iteration],
806 localDiffusionsGrid.
L, localDiffusionsGrid.
pc, localDiffusionsGrid.
alpha, localDiffusionsGrid.
Jacobian,
807 pc_LowerBoundaryCondition.
arr,
808 pc_UpperBoundaryCondition.
arr,
809 alpha_LowerBoundaryCondition.
arr,
810 alpha_UpperBoundaryCondition.
arr,
821 if (localDiffusionsGrid.
pc.
size > 1 && localDiffusionsGrid.
alpha.
size > 1) {
822 Output::echo(3,
"Calculating local diffusions using Kyung-Chan's approximation... \n");
825 parameters.
Lpp[iteration],
829 localDiffusionsGrid.
L, localDiffusionsGrid.
pc, localDiffusionsGrid.
alpha, localDiffusionsGrid.
Jacobian,
830 pc_LowerBoundaryCondition.
arr,
831 pc_UpperBoundaryCondition.
arr,
832 alpha_LowerBoundaryCondition.
arr,
833 alpha_UpperBoundaryCondition.
arr,
848 throw error_msg(
"APPROX_ERROR",
"Approximation unknown");
865 SL.
Update(iteration, parameters.
SL, localDiffusionsGrid);
869 localDiffusionsGrid.
L, localDiffusionsGrid.
pc, localDiffusionsGrid.
alpha,
872 parameters.
Lpp[iteration],
873 parameters.
tau[iteration],
874 parameters.
tauLpp[iteration], parameters.
Kp[iteration]);
890 for (il = 0; il < localDiffusionsGrid.
L.
size; il++)
891 for (im = 0; im < localDiffusionsGrid.
pc.
size; im++)
892 for (ia = 0; ia < localDiffusionsGrid.
alpha.
size; ia++)
911 double LaaLR, LaaRL, LmmLR, LmmRL;
912 double LmaLL, LmaLR, LmaRL, LmaRR;
913 double LamLL, LamLR, LamRL, LamRR;
917 double maxerror = 0, aveerror = 0;
919 for (il = 0; il < localDiffusionsGrid.
L.
size; il++) {
920 for (im = 0; im < localDiffusionsGrid.
pc.
size; im++) {
922 PSD_Der_a_R[il][im][ia] = localDiffusionsGrid.
Jacobian[il][im][ia] * Dpca.
CurrentDxx.
arr[il][im][ia] * (psdLocalDiffusions.
arr[il][im][ia+1] - psdLocalDiffusions.
arr[il][im][ia])/(localDiffusionsGrid.
alpha.
arr[il][im][ia+1] - localDiffusionsGrid.
alpha.
arr[il][im][ia]);
923 for (ia = 1; ia < localDiffusionsGrid.
alpha.
size-1; ia++) {
924 PSD_Der_a_L[il][im][ia] = localDiffusionsGrid.
Jacobian[il][im][ia] * Dpca.
CurrentDxx.
arr[il][im][ia] * (psdLocalDiffusions.
arr[il][im][ia] - psdLocalDiffusions.
arr[il][im][ia-1])/(localDiffusionsGrid.
alpha.
arr[il][im][ia] - localDiffusionsGrid.
alpha.
arr[il][im][ia-1]);
925 PSD_Der_a_R[il][im][ia] = localDiffusionsGrid.
Jacobian[il][im][ia] * Dpca.
CurrentDxx.
arr[il][im][ia] * (psdLocalDiffusions.
arr[il][im][ia+1] - psdLocalDiffusions.
arr[il][im][ia])/(localDiffusionsGrid.
alpha.
arr[il][im][ia+1] - localDiffusionsGrid.
alpha.
arr[il][im][ia]);
928 PSD_Der_a_L[il][im][ia] = localDiffusionsGrid.
Jacobian[il][im][ia] * Dpca.
CurrentDxx.
arr[il][im][ia] * (psdLocalDiffusions.
arr[il][im][ia] - psdLocalDiffusions.
arr[il][im][ia-1])/(localDiffusionsGrid.
alpha.
arr[il][im][ia] - localDiffusionsGrid.
alpha.
arr[il][im][ia-1]);
930 for (ia = 0; ia < localDiffusionsGrid.
alpha.
size; ia++) {
932 PSD_Der_m_R[il][im][ia] = localDiffusionsGrid.
Jacobian[il][im][ia] * Dpca.
CurrentDxx.
arr[il][im][ia]*(psdLocalDiffusions.
arr[il][im+1][ia] - psdLocalDiffusions.
arr[il][im][ia])/(localDiffusionsGrid.
pc.
arr[il][im+1][ia] - localDiffusionsGrid.
pc.
arr[il][im][ia]);
933 for (im = 1; im < localDiffusionsGrid.
pc.
size-1; im++) {
934 PSD_Der_m_L[il][im][ia] = localDiffusionsGrid.
Jacobian[il][im][ia] * Dpca.
CurrentDxx.
arr[il][im][ia]*(psdLocalDiffusions.
arr[il][im][ia] - psdLocalDiffusions.
arr[il][im-1][ia])/(localDiffusionsGrid.
pc.
arr[il][im][ia] - localDiffusionsGrid.
pc.
arr[il][im-1][ia]);
935 PSD_Der_m_R[il][im][ia] = localDiffusionsGrid.
Jacobian[il][im][ia] * Dpca.
CurrentDxx.
arr[il][im][ia]*(psdLocalDiffusions.
arr[il][im+1][ia] - psdLocalDiffusions.
arr[il][im][ia])/(localDiffusionsGrid.
pc.
arr[il][im+1][ia] - localDiffusionsGrid.
pc.
arr[il][im][ia]);
937 im = localDiffusionsGrid.
pc.
size-1;
938 PSD_Der_m_L[il][im][ia] = localDiffusionsGrid.
Jacobian[il][im][ia] * Dpca.
CurrentDxx.
arr[il][im][ia]*(psdLocalDiffusions.
arr[il][im][ia] - psdLocalDiffusions.
arr[il][im-1][ia])/(localDiffusionsGrid.
pc.
arr[il][im][ia] - localDiffusionsGrid.
pc.
arr[il][im-1][ia]);
942 for (il = 0; il < localDiffusionsGrid.
L.
size; il++) {
943 for (im = 1; im < localDiffusionsGrid.
pc.
size-1; im++) {
944 for (ia = 1; ia < localDiffusionsGrid.
alpha.
size-1; ia++) {
947 LaaLR = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
948 (PSD_Der_a_R[il][im][ia] - PSD_Der_a_R[il][im][ia-1])/
949 (localDiffusionsGrid.
alpha.
arr[il][im][ia] - localDiffusionsGrid.
alpha.
arr[il][im][ia-1]);
950 LaaRL = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
951 (PSD_Der_a_L[il][im][ia+1] - PSD_Der_a_L[il][im][ia])/
952 (localDiffusionsGrid.
alpha.
arr[il][im][ia+1] - localDiffusionsGrid.
alpha.
arr[il][im][ia]);
955 LmmLR = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
956 (PSD_Der_m_R[il][im][ia] - PSD_Der_m_R[il][im-1][ia])/
957 (localDiffusionsGrid.
pc.
arr[il][im][ia] - localDiffusionsGrid.
pc.
arr[il][im-1][ia]);
958 LmmRL = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
959 (PSD_Der_m_L[il][im+1][ia] - PSD_Der_m_L[il][im][ia])/
960 (localDiffusionsGrid.
pc.
arr[il][im+1][ia] - localDiffusionsGrid.
pc.
arr[il][im][ia]);
963 LmaRR = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
964 (PSD_Der_a_R[il][im+1][ia] - PSD_Der_a_R[il][im][ia])/
965 (localDiffusionsGrid.
pc.
arr[il][im+1][ia] - localDiffusionsGrid.
pc.
arr[il][im][ia]);
966 LmaRL = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
967 (PSD_Der_a_L[il][im+1][ia] - PSD_Der_a_L[il][im][ia])/
968 (localDiffusionsGrid.
pc.
arr[il][im+1][ia] - localDiffusionsGrid.
pc.
arr[il][im][ia]);
969 LmaLR = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
970 (PSD_Der_a_R[il][im][ia] - PSD_Der_a_R[il][im-1][ia])/
971 (localDiffusionsGrid.
pc.
arr[il][im][ia] - localDiffusionsGrid.
pc.
arr[il][im-1][ia]);
972 LmaLL = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
973 (PSD_Der_a_L[il][im][ia] - PSD_Der_a_L[il][im-1][ia])/
974 (localDiffusionsGrid.
pc.
arr[il][im][ia] - localDiffusionsGrid.
pc.
arr[il][im-1][ia]);
978 LamRR = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
979 (PSD_Der_m_R[il][im][ia+1] - PSD_Der_m_R[il][im][ia])/
980 (localDiffusionsGrid.
alpha.
arr[il][im][ia+1] - localDiffusionsGrid.
alpha.
arr[il][im][ia]);
981 LamRL = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
982 (PSD_Der_m_L[il][im][ia+1] - PSD_Der_m_L[il][im][ia])/
983 (localDiffusionsGrid.
alpha.
arr[il][im][ia+1] - localDiffusionsGrid.
alpha.
arr[il][im][ia]);
984 LamLR = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
985 (PSD_Der_m_R[il][im][ia] - PSD_Der_m_R[il][im][ia-1])/
986 (localDiffusionsGrid.
alpha.
arr[il][im][ia] - localDiffusionsGrid.
alpha.
arr[il][im][ia-1]);
987 LamLL = 1.0/(localDiffusionsGrid.
Jacobian[il][im][ia]) *
988 (PSD_Der_m_L[il][im][ia] - PSD_Der_m_L[il][im][ia-1])/
989 (localDiffusionsGrid.
alpha.
arr[il][im][ia] - localDiffusionsGrid.
alpha.
arr[il][im][ia-1]);
991 error[il][im][ia] = psdLocalDiffusions.
arr[il][im][ia] - psdPrevious[il][im][ia]
992 - (LaaLR + LaaRL) / 2 * parameters.
timeStep
993 - (LmmLR + LmmRL) / 2 * parameters.
timeStep
994 - (LmaRR + LmaRL + LmaLR + LmaLL) / 4 * parameters.
timeStep
995 - (LamRR + LamRL + LamLR + LamLL) / 4 * parameters.
timeStep;
999 error[il][im][ia] += psdLocalDiffusions.
arr[il][im][ia]/taulc * parameters.
timeStep;
1002 aveerror += fabs(error[il][im][ia]);
1003 maxerror = (fabs(maxerror) > fabs(error[il][im][ia]))?maxerror:error[il][im][ia];
1007 aveerror = aveerror / (localDiffusionsGrid.
L.
size) / (localDiffusionsGrid.
epc.
size-2) / (localDiffusionsGrid.
alpha.
size-2);
1008 Output::echo(20,
"Max calculation error = %e; ave calculation error = %e\n", maxerror, aveerror);
1010 stringstream filename;
1011 filename <<
"./Debug/Calc_error_" << iteration <<
".plt";
1015 psdPrevious = psdLocalDiffusions.
arr;
1026 Output::echo(0,
"Step %d of %d. Time = %d days, %.0f hours. ", iteration, parameters.
totalIterationsNumber,
int(simulation_time), (simulation_time-
int(simulation_time))*24);
1030 Output1DValues(output1D, Daa, Dpcpc, Dpca, DaaLpp, DpcpcLpp, DpcaLpp, simulation_time, parameters, iteration);
1052 if (radialDiffusionGrid.
L.
size > 1) {
1066 psdRadialDiffusion.
Interpolate(psdLocalDiffusions, parameters.
interpolation, localDiffusionsGrid, radialDiffusionGrid, pc_LowerBoundaryCondition.
arr, psdLocalDiffusions.
arr.
ySlice(radialDiffusionGrid.
pc.
size-1));
1078 current_time = std::clock();
1080 Output::echo(0,
"Total time %0.2f sec\n", ( current_time - zero_time ) / (
double)CLOCKS_PER_SEC );
1097 }
catch (exception
exp) {
1098 Output::echo(0,
"\nCRITICAL ERROR:\n%s\n", exp.what());
Holds list of instances of DiffusionCoefficient class of same type (like Daa, Dpp, etc), but produced by different waves (Daa_chorus, Daa_EMIC, etc).
void Update(int iteration, Parameters_structure::SL_structure &SL, Grid &grid)
GridElement alpha
grid elements to be used
Matrix3D< double > arr
array of grid points
Matrix1D< double > Kp
Kp array.
void Diffusion_L(double dt, double Lpp, DiffusionCoefficient &DLL, GridElement &L, GridElement &pc, GridElement &alpha, Matrix3D< double > Jacobian, Matrix2D< double > lowerBoundaryCondition, Matrix2D< double > upperBoundaryCondition, string lowerBoundaryCondition_calculationType, string upperBoundaryCondition_calculationType, double tau, double tauLpp)
void Interpolate(PSD &otherPSD, Parameters_structure::Interpolation interpolationParameters, Grid &oldGrid, Grid &newGrid, Matrix2D< double > newGrid_pc_lowerBoundaryCondition, Matrix2D< double > newGrid_pc_upperBoundaryCondition)
Interpolation function.
string folderName
output folder name.
struct Parameters_structure::GridElement radialDiffusionGrid_L
grid elements for local and radial diffusion, and corresponding types of elements ...
bool ActivateAndScale(double time, double Kp)
void SourcesAndLosses(Parameters_structure parameters, GridElement &L, GridElement &pc, GridElement &alpha, AdditionalSourcesAndLosses &SL, double dt, double Lpp, double tau, double tauLpp, double Kp)
Calculates any additional sources and losses due to magnetopause.
struct Parameters_structure::Interpolation interpolation
Interpolation parameters.
Matrix2D< T > ySlice(int p_y) const
bool useLossCone
specify whether using loss cone or not
#define __VERB_VERSION_NUMBER__
Number of current VERB code version.
void Output1DValues(ofstream &output1D, DiffusionCoefficientsGroup &Daa, DiffusionCoefficientsGroup &Dpcpc, DiffusionCoefficientsGroup &Dpca, DiffusionCoefficientsGroup &DaaLpp, DiffusionCoefficientsGroup &DpcpcLpp, DiffusionCoefficientsGroup &DpcaLpp, double &time, Parameters_structure ¶meters, int iteration)
Function used in main() to print values of diffusion coefficients in output file. ...
void Diffusion_alpha(double dt, double Lpp, DiffusionCoefficient &Daa, DiffusionCoefficient &DaaLpp, GridElement &L, GridElement &pc, GridElement &alpha, Matrix3D< double > Jacobian, Matrix2D< double > alpha_lowerBoundaryCondition, Matrix2D< double > alpha_upperBoundaryCondition, string alpha_lowerBoundaryCondition_calculationType, string alpha_upperBoundaryCondition_calculationType)
void writeToFile(string filename)
Save matrix to a file.
Holds diffusion coefficient matrix and routines to load and calculate it.
bool useAlphaDiffusion
Using diffusions flags.
void Diffusion_pc_alpha(double dt, double Lpp, DiffusionCoefficient &Dpcpc, DiffusionCoefficient &DpcpcLpp, DiffusionCoefficient &Daa, DiffusionCoefficient &DaaLpp, DiffusionCoefficient &Dpca, DiffusionCoefficient &DpcaLpp, GridElement &L, GridElement &pc, GridElement &alpha, Matrix3D< double > Jacobian, Matrix2D< double > pc_lowerBoundaryCondition, Matrix2D< double > pc_upperBoundaryCondition, Matrix2D< double > alpha_lowerBoundaryCondition, Matrix2D< double > alpha_upperBoundaryCondition, string pc_lowerBoundaryCondition_calculationType, string pc_upperBoundaryCondition_calculationType, string alpha_lowerBoundaryCondition_calculationType, string alpha_upperBoundaryCondition_calculationType)
string calculationType
Type: constant function/constant derivative.
void Create_Initial_PSD(Parameters_structure::PSD parameters, Grid &grid, BoundaryCondition L_UpperBoundaryCondition)
Create inital PSD (Steady State)
void Output_without_grid(double time)
void set_output_lvl(int new_outputLvl)
Set the verbose level until redefine.
struct Parameters_structure::General_Output_parameters general_Output_parameters
instance of the struct
Matrix2D< double > arr
2D matrix of boundary conditions
DiffusionCoefficient CurrentDxx
flag, indicated that the initialization was passed
Matrix1D< double > Bf
Boundary flux array.
static const double exp
exponential
void echo(int msglvl, const char *format,...)
Basic function! Write the message to the file.
void open_log_file(string filename)
Open file stream for log file.
void AllocateMemory(int size_x, int size_y, int size_z)
void flush_log_file()
flush file stream for log file
string radialDiffusionGrid_filename
filename, if grid load from file
struct Parameters_structure::SL_structure SL
instance of the struct
void Create_Grid(Parameters_structure::GridElement parameters_L, Parameters_structure::GridElement parameters_pc, Parameters_structure::GridElement parameters_alpha, Parameters_structure::GridElement parameters_epc, string grid_filename, string gridType, Grid SecondGrid=Grid())
int main(int argc, char *argv[])
Main code.
int iterStep
Number of iterations between outputs.
string approximationMethod
Approximation method. Check StrToVal(string input, ApproximationMethods &place) for known values...
void MakeBoundaryCondition(Parameters_structure::BoundaryCondition parameters, Matrix2D< double > psd2DSlice, Matrix2D< double > gridElement1, Matrix2D< double > gridElement2)
void DiffusionMixTermExplicit(double dt, double Lpp, DiffusionCoefficient &Dpca, DiffusionCoefficient &DpcaLpp, GridElement &L, GridElement &pc, GridElement &alpha, Matrix3D< double > Jacobian, Matrix2D< double > pc_lowerBoundaryCondition, Matrix2D< double > pc_upperBoundaryCondition, Matrix2D< double > alpha_lowerBoundaryCondition, Matrix2D< double > alpha_upperBoundaryCondition, string pc_lowerBoundaryCondition_calculationType, string pc_upperBoundaryCondition_calculationType, string alpha_lowerBoundaryCondition_calculationType, string alpha_upperBoundaryCondition_calculationType)
Matrix1D< double > Lpp
Lpp array.
string DLLType
DLL type (which method to use to calculate). Check StrToVal(string input, DLLTypes &place) for known ...
string type
Type. Check StrToVal(string input, BoundaryConditionTypes &place) for known values.
void Diffusion_pc(double dt, double Lpp, DiffusionCoefficient &Dpcpc, DiffusionCoefficient &DpcpcLpp, GridElement &L, GridElement &pc, GridElement &alpha, Matrix3D< double > Jacobian, Matrix2D< double > pc_lowerBoundaryCondition, Matrix2D< double > pc_upperBoundaryCondition, string pc_lowerBoundaryCondition_calculationType, string pc_upperBoundaryCondition_calculationType)
Main parameters structure that holds smaller structures for individual parameters.
bool useEnergyDiffusion
Using diffusions flags.
int size
size of grid element
bool useEnergyAlphaMixedTerms
int totalIterationsNumber
Total number of iterations.
Matrix3D< double > Jacobian
void close_log_file()
Close file stream for log file.
double bounce_time_new(double L, double pc, double alpha)
bounce time
void Create_SL(Parameters_structure::SL_structure &SL, Grid &grid, int numberOfIterations, double timeStep)
Matrix3D< double > arr
array of diffusion coefficients
bool NoNegative
Artificially no negative PSD.
bool Load_parameters(string filename)
void CreateAllDiffusionCoefficients(DiffusionCoefficient &DLL, DiffusionCoefficientsGroup &Daa, DiffusionCoefficientsGroup &Dpcpc, DiffusionCoefficientsGroup &Dpca, DiffusionCoefficientsGroup &DaaLpp, DiffusionCoefficientsGroup &DpcpcLpp, DiffusionCoefficientsGroup &DpcaLpp, Parameters_structure ¶meters, Grid &radialDiffusionGrid, Grid &localDiffusionsGrid)
string radialDiffusionGrid_type
Local diffusions grid type.
void Output1DHeaders(ofstream &output1D, DiffusionCoefficientsGroup &Daa, DiffusionCoefficientsGroup &Dpcpc, DiffusionCoefficientsGroup &Dpca, DiffusionCoefficientsGroup &DaaLpp, DiffusionCoefficientsGroup &DpcpcLpp, DiffusionCoefficientsGroup &DpcaLpp)
Function used in main() to print header in output file. Header needs to understand the output file...
void Diffusion_pc_alpha_KC(double dt, double Lpp, DiffusionCoefficient &Dpcpc, DiffusionCoefficient &DpcpcLpp, DiffusionCoefficient &Daa, DiffusionCoefficient &DaaLpp, DiffusionCoefficient &Dpca, DiffusionCoefficient &DpcaLpp, GridElement &L, GridElement &pc, GridElement &alpha, Matrix3D< double > Jacobian, Matrix2D< double > pc_lowerBoundaryCondition, Matrix2D< double > pc_upperBoundaryCondition, Matrix2D< double > alpha_lowerBoundaryCondition, Matrix2D< double > alpha_upperBoundaryCondition, string pc_lowerBoundaryCondition_calculationType, string pc_upperBoundaryCondition_calculationType, string alpha_lowerBoundaryCondition_calculationType, string alpha_upperBoundaryCondition_calculationType)
bool outputModelMatrix
Output model matrix.
Matrix3D< double > arr
array of PSD values
Phase Space Density class.
GridElement L
grid elements to be used
double timeStep
Time step, hours.
double alc(double L)
Loss cone calculations.
Holds upper and lower boundary conditions.
Matrix2D< T > xSlice(int p_x) const
Parameters_structure parameters
Parameters structure, with all parameters from the parameters.ini file. The default parameters define...
GridElement epc
grid elements to be used
Error message - stack of single_errors.
string localDiffusionsGrid_filename
filename, if grid load from file
Matrix1D< double > tauLpp
Lifetime inside of the plasmasphere.
void Output(string filename)
bool check_time(int iter, int d_iter)
Small routine for checking if it is right iteration for the next output.
string fileName1D
1d-output file name.
string localDiffusionsGrid_type
Radial diffusions grid type.
Matrix1D< double > tau
Lifetime out of the plasmasphere.
Computational grid composed of 3 different GridElement.
void writeToFile(string filename)
static const double zero_f
Minimum value for PSD.
CalculationMatrix matr_C
calculation matrices A, B, and C
CalculationMatrix matr_B
calculation matrices A, B, and C
void Update(int iteration, Matrix2D< double > PSD_2D_Slice, double time=-1, double dt=-1)
bool useRadialDiffusion
Using diffusions flags.
void LoadBoundaryCondition(Parameters_structure::BoundaryCondition parameters, Matrix2D< double > gridElement1, Matrix2D< double > gridElement2)
struct Parameters_structure::PSD psdRadialDiffusion
PSD parameters for radial diffusion.
int outputLvl
Detalization level of screen output.
void MakeDLL(double Kp)
Making DLL.
Matrix2D< T > zSlice(int p_z) const
struct Parameters_structure::PSD psdLocalDiffusions
PSD parameters for local diffusions.
GridElement pc
grid elements to be used
CalculationMatrix matr_A
calculation matrices A, B, and C