Commit 8b25202b authored by niklas.baumgarten's avatar niklas.baumgarten
Browse files

added coarse level, added flux plot added time sereies

parent fd02b33a
...@@ -18,29 +18,30 @@ Problem = StochasticPollution2D ...@@ -18,29 +18,30 @@ Problem = StochasticPollution2D
StochasticField = LogNormal StochasticField = LogNormal
#StochasticField = Gaussian #StochasticField = Gaussian
#Experiment = ConvergenceTest Experiment = ConvergenceTest
Experiment = MLMCExperiment #Experiment = MLMCExperiment
#Experiment = MLMCOverEpsilon #Experiment = MLMCOverEpsilon
initLevels = 4, 5, initLevels = 4, 5
initSampleAmount = 1, 12, initSampleAmount = 0, 1
#l_init = 3, 4, 5 #l_init = 3, 4, 5
#M_init = 45, 15, 5 #M_init = 45, 15, 5
epsilon = 0.01 epsilon = 0.003
mcOnly = false mcOnly = false
epsilon_lst = 0.05, 0.01, 0.005 epsilon_lst = 0.05, 0.01, 0.005
maxLevel = 5 maxLevel = 8
plevel = 2 plevel = 2
degree = 2 degree = 2
functional = L2 #functional = L2
#functional = Energy functional = Energy
#funcitonal = H1 #funcitonal = H1
#functional = Outflow #functional = Outflow
#functional = Goal #functional = Goal
#functional = Mass
MCVerbose = 2 MCVerbose = 2
MLMCVerbose = 2 MLMCVerbose = 2
...@@ -80,11 +81,12 @@ MuteLevel = -1 ...@@ -80,11 +81,12 @@ MuteLevel = -1
NewtonVerbose = -1 NewtonVerbose = -1
DebugLevel = -1 DebugLevel = -1
#Transport #Transport
startTime = 0.0
endTime = 2.0
TimeSeries = uniform TimeSeries = uniform
plot_tStep = 1 T = 2
T = 1.6
dt = 0.005 dt = 0.005
Kmax = 250 Kmax = 250
Keps = 1e-5 Keps = 1e-5
...@@ -96,4 +98,7 @@ Keps = 1e-5 ...@@ -96,4 +98,7 @@ Keps = 1e-5
rkorder = -2 #Implicit MP-rule rkorder = -2 #Implicit MP-rule
#rkorder = 1000 #Polynomial Krylov-Arnoldi method #rkorder = 1000 #Polynomial Krylov-Arnoldi method
gamma = 0.01 gamma = 0.01
\ No newline at end of file
plotPermForTransport = true
\ No newline at end of file
...@@ -19,7 +19,7 @@ public: ...@@ -19,7 +19,7 @@ public:
const char *Name() const override { return "DGTransportAssemble"; } const char *Name() const override { return "DGTransportAssemble"; }
virtual double Residual(const Vector &u, Vector &b) const { double Residual(const Vector &u, Vector &b) const override {
b = 0; b = 0;
Matrix *massMatrix = new Matrix(u); Matrix *massMatrix = new Matrix(u);
...@@ -302,7 +302,7 @@ public: ...@@ -302,7 +302,7 @@ public:
energy += w * (U * U); energy += w * (U * U);
} }
} }
return 0.5 * PPM->Sum(energy); //TODO was soll das sein? M.E. L2 ! return 0.5 * PPM->Sum(energy);
} }
double Error(double t, const Vector &u) const { double Error(double t, const Vector &u) const {
......
...@@ -38,12 +38,19 @@ void MonteCarloTransport::Method() { ...@@ -38,12 +38,19 @@ void MonteCarloTransport::Method() {
plot.vtk_cellvector(name); plot.vtk_cellvector(name);
} }
solveTransport(l, fineQ, fineCost, fineSolution); solveTransport(l, true, fineQ, fineCost, fineSolution);
if (!onlyFineLevel) { if (!onlyFineLevel) {
stochasticField->GetCoarseSample(l - 1, fineFlux, coarseFlux); stochasticField->GetCoarseSample(l - 1, fineFlux, coarseFlux);
assemble->problem->LoadNewSample(&coarseFlux); assemble->problem->LoadNewSample(&coarseFlux);
solveTransport(l, fineQ, fineCost, fineSolution);
if (plotting > 1) {
const char *name = buildString(m, "flux_down");
plot.celldata(coarseFlux, 3);
plot.vtk_cellvector(name);
}
solveTransport(l, false, coarseQ, coarseCost, coarseSolution);
} else { } else {
coarseQ = 0.0, coarseCost = 0.0, coarseSolutionUp = 0.0; coarseQ = 0.0, coarseCost = 0.0, coarseSolutionUp = 0.0;
} }
...@@ -57,17 +64,37 @@ void MonteCarloTransport::Method() { ...@@ -57,17 +64,37 @@ void MonteCarloTransport::Method() {
logger->endMethodMSG(); logger->endMethodMSG();
} }
void MonteCarloTransport::solveTransport(int l, double &valueQ, double &cost, void MonteCarloTransport::solveTransport(int l,
bool fineLevel,
double &valueQ,
double &cost,
Vector &solution) { Vector &solution) {
Date startSolving; // Date startSolving;
Solver solver; Solver solver;
Preconditioner *BJ = GetPC("PointBlockJacobi"); Preconditioner *BJ = GetPC("PointBlockJacobi");
Solver IM(BJ); Solver IM(BJ);
DGTimeIntegrator timeIntegrator(solver, solver, solver, IM); DGTimeIntegrator timeIntegrator(solver, solver, solver, IM);
assemble->PrintInfo(solution); TimeSeries timeSeries = getTimeSeries(fineLevel);
assemble->resetTAssemble();
timeIntegrator(*assemble, timeSeries, solution); timeIntegrator(*assemble, timeSeries, solution);
cost = solution.pSize() * assemble->getStep();
if (functional == "Mass") valueQ = assemble->Mass(solution);
if (functional == "Energy") valueQ = assemble->Energy(solution);
if (functional == "Outflow")
valueQ = assemble->InFlowOutFlowRate(solution).second;
logger->logMSGV2("Solving on level " + to_string(l) + " took " + // logger->logMSGV2("Solving on level " + to_string(l) + " took " +
to_string((startSolving - Date()).t) + " seconds"); // to_string((startSolving - Date()).t) + " seconds");
}
TimeSeries MonteCarloTransport::getTimeSeries(bool fineLevel) {
int scaling = 8;
if (fineLevel)
return TimeSeries(startTime, endTime,
scaling * (int) (pow(2, l) * (endTime - startTime)));
else
return TimeSeries(startTime, endTime,
scaling * (int) (pow(2, l - 1) * (endTime - startTime)));
} }
...@@ -3,23 +3,29 @@ ...@@ -3,23 +3,29 @@
#include "MonteCarlo.h" #include "MonteCarlo.h"
#include "assemble/DGTransportAssemble.h" #include "assemble/DGTransportAssemble.h"
#include "TimeIntegrator.h" #include "timestepping/TimeIntegrator.h"
class MonteCarloTransport : public MonteCarlo { class MonteCarloTransport : public MonteCarlo {
private: private:
Plot plot; Plot plot;
void solveTransport(int l, double &valueQ, double &cost, Vector &solution); void solveTransport(int l,
bool fineLevel,
double &valueQ,
double &cost,
Vector &solution);
TimeSeries getTimeSeries(bool fineLevel);
public: public:
DGTransportAssemble *assemble; DGTransportAssemble *assemble;
TimeSeries timeSeries; double startTime = 0.0;
double endTime = 0.0;
MatrixGraphs fluxMatrixGraphs; MatrixGraphs fluxMatrixGraphs;
CellMatrixGraphs solMatrixGraphs; CellMatrixGraphs solMatrixGraphs;
Transfer *transfer; Transfer *transfer;
MonteCarloTransport(int l, int dM, bool baseLevel, Meshes *meshes, MonteCarloTransport(int l, int dM, bool baseLevel, Meshes *meshes,
...@@ -30,6 +36,8 @@ public: ...@@ -30,6 +36,8 @@ public:
solMatrixGraphs(*meshes, *assemble->disc), solMatrixGraphs(*meshes, *assemble->disc),
plot(Plot((*meshes)[l], meshes->dim())) { plot(Plot((*meshes)[l], meshes->dim())) {
config.get("startTime", startTime);
config.get("endTime", endTime);
// transfer = new MatrixTransfer(*assemble); // transfer = new MatrixTransfer(*assemble);
// transfer->Construct((solMatrixGraphs)[l - pLevel], (solMatrixGraphs)[l - pLevel - 1]); // transfer->Construct((solMatrixGraphs)[l - pLevel], (solMatrixGraphs)[l - pLevel - 1]);
} }
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment