Commit be4251e0 authored by niklas.baumgarten's avatar niklas.baumgarten
Browse files

Merge branch '36-further-sc-tests' into 'feature'

Resolve "Further SC Tests"

Closes #36

See merge request !55
parents 1bbb666e 5117adea
Pipeline #168478 passed with stages
in 18 minutes and 52 seconds
......@@ -29,10 +29,14 @@ void StochasticCollocation::method() {
while (aggregate.ctr.dMcomm != 0) {
updateSampleIds(fineSolution, coarseSolution);
// Todo implement coarse solution solving
pdeSolver->DrawSample(fineId);
pdeSolver->Run(fineSolution);
aggregate.Update(fineSolution, coarseSolution);
// Todo: implement error estimation!
// errors.Estimate(aggregate);
}
}
......
......@@ -5,6 +5,9 @@
#include <utility>
// Todo
// * this class should be merged with MC into Estimator class.
// * Then provide purely abstract Estimator Interface
class StochasticCollocation : public Estimator {
protected:
......
......@@ -17,6 +17,8 @@ protected:
int outputs;
//double alpha;
double weight = 0.0;
RVector sample{};
......@@ -37,6 +39,23 @@ protected:
};
public:
// Todo:
// * Use other grid types. (e.g. for normally distributed rv)
// * If so, write tests and split this class
// * Maybe overload other generators than SampleGenerator<RVector>
//
// Todo:
// * Find out which level corresponds to which total degree and exact degree
// * Which depth_rule will be the most useful for our purpose
//
// Todo:
// * Move functions to cpp and only include tasmanian there!
//
// Todo:
// * Distribute samples ob different processes
SparseGridGenerator(const Meshes &meshes, int dimension, int outputs, int level,
TasGrid::TypeDepth depth = TasGrid::type_level,
TasGrid::TypeOneDRule rule = TasGrid::rule_clenshawcurtis) :
......@@ -59,6 +78,9 @@ public:
return weights;
}
// Todo implement linear trafo for areas other than [-1, 1]^N
// Todo Use TASMANIAN interface for that
double SampleWeight(const SampleID &id) {
return weights[id.number];
}
......@@ -68,6 +90,7 @@ public:
}
double SumOfWeights() const {
// Todo has to be adapted to area
return 2.0 * dimension;
}
......
......@@ -23,5 +23,15 @@ CreateStochasticDummyProblem(const std::string &problemName, const Meshes &meshe
return new SparseGrid1DGeneratorProblem(meshes);
if (problemName == "SparseGrid2DGeneratorProblem")
return new SparseGrid2DGeneratorProblem(meshes);
if (problemName == "SparseGrid2DPolyDeg2Problem")
return new SparseGrid2DPolyDeg2Problem(meshes);
if (problemName == "SparseGrid2DPD2qptotProblem")
return new SparseGrid2DPD2qptotProblem(meshes);
if (problemName == "SparseGrid2DHermiteProblem")
return new SparseGrid2DHermiteProblem(meshes);
if (problemName == "SparseGrid2DPD2qptotProblemFail")
return new SparseGrid2DPD2qptotProblemFail(meshes);
if (problemName == "SparseGrid2DPD3qptotProblemFail")
return new SparseGrid2DPD3qptotProblemFail(meshes);
Exit(problemName + " not found")
}
\ No newline at end of file
......@@ -8,7 +8,7 @@
#include "SparseGridGenerator.hpp"
#include "HybridFluxGenerator.hpp"
#include "CirculantEmbedding.hpp"
#include "math.h"
class IStochasticProblem {
protected:
......@@ -249,6 +249,96 @@ public:
}
};
class SparseGrid2DPolyDeg2Problem : public SparseGridGeneratorProblem {
public:
explicit SparseGrid2DPolyDeg2Problem(const Meshes &meshes) :
SparseGridGeneratorProblem(meshes, SparseGridGenerator(meshes, 2, 0, 4)) {}
double FunctionEvaluation() override {
RVector sample = this->generator.EvalSample();
return 3.0 * (1 + sample[0] + pow(sample[0],2) + sample[1] + pow(sample[1],2) + sample[0] * sample[1]);
}
string Name() const override {
return "SparseGrid2DPolyDeg2Problem";
}
};
//integrating polynomials upto given exact degree n or less
class SparseGrid2DPD2qptotProblem : public SparseGridGeneratorProblem {
public:
explicit SparseGrid2DPD2qptotProblem(const Meshes &meshes) :
SparseGridGeneratorProblem(meshes, SparseGridGenerator(meshes, 2, 0, 2, TasGrid::type_qptotal)) {}
double FunctionEvaluation() override {
RVector sample = this->generator.EvalSample();
double y1 = sample[0];
double y2 = sample[1];
return 3.0 * (1 + y1 + pow(y1,2) + y2 + pow(y1,2) + y1 * y2);
}
string Name() const override {
return "SparseGrid2DPD2qptotProblem";
}
};
class SparseGrid2DHermiteProblem : public SparseGridGeneratorProblem {
public:
explicit SparseGrid2DHermiteProblem(const Meshes &meshes) :
SparseGridGeneratorProblem(meshes, SparseGridGenerator(meshes, 1, 0, 8,
TasGrid::type_qptotal, TasGrid::rule_gausshermite)) {}
double FunctionEvaluation() override {
RVector sample = this->generator.EvalSample();
//double y1 = sample[0];
return 1 / sqrt(M_PI);
}
string Name() const override {
return "SparseGrid2DhermiteProblem";
}
};
// this test is meant to fail. We are using a sparse grid that only integrates functions exactly upto degree 1 or less
class SparseGrid2DPD2qptotProblemFail : public SparseGridGeneratorProblem {
public:
explicit SparseGrid2DPD2qptotProblemFail(const Meshes &meshes) :
SparseGridGeneratorProblem(meshes, SparseGridGenerator(meshes, 2, 0, 2, TasGrid::type_qptotal)) {}
double FunctionEvaluation() override {
RVector sample = this->generator.EvalSample();
double y1 = sample[0];
double y2 = sample[1];
return 1 + y1 + pow(y1,4) + y2 + pow(y2,5) + y1 * y2;
}
string Name() const override {
return "SparseGrid2DPD2qptotProblemFail";
}
};
class SparseGrid2DPD3qptotProblemFail : public SparseGridGeneratorProblem {
public:
explicit SparseGrid2DPD3qptotProblemFail(const Meshes &meshes) :
SparseGridGeneratorProblem(meshes, SparseGridGenerator(meshes, 2, 0, 4, TasGrid::type_qptotal)) {}
double FunctionEvaluation() override {
RVector sample = this->generator.EvalSample();
double y1 = sample[0];
double y2 = sample[1];
return pow(y1,4) + pow(y2,6);
}
string Name() const override {
return "SparseGrid2DPD3qptotProblemFail";
}
};
// Todo find more test cases in tasmanian examples and write sparse grid problems
StochasticDummyProblem *
CreateStochasticDummyProblem(const std::string &problemName, const Meshes &meshes);
......
......@@ -154,5 +154,7 @@ CreateStochasticEllipticProblem(const string &problemName, const Meshes &meshes)
return new Laplace2D(meshes);
if (problemName == "StochasticLaplace2DTest")
return new StochasticLaplace2DTest(meshes);
if (problemName == "SparseGridTestPDE")
return new SparseGridTestPDE(meshes);
Exit(problemName + " not found")
}
\ No newline at end of file
......@@ -23,6 +23,40 @@ public:
}
};
class SparseGridTestPDE : public IStochasticEllipticProblem {
UniformDistributionRVector rVectorGenerator;
private:
double a_min = 1/100 ;
public:
explicit SparseGridTestPDE(const Meshes &meshes) : IStochasticEllipticProblem(meshes),
rVectorGenerator(UniformDistributionRVector(meshes, 4, 0.0, 1.0)) {}
void DrawSample(const SampleID &id) override {
rVectorGenerator.DrawSample(id);
}
//Sets dirichlet boundary condition
Scalar Solution(const Point &x) const override{
return 0.0;
}
//sets neumann boundary conditions
VectorField Flux(const Point &x) const override {
return VectorField(1.0,0.0);
}
Scalar Load(const Point &x) const override {
return 0.0;
}
Tensor Permeability(const cell &c) const override {
return a_min + exp((cos(M_PI * c()[1]) * this -> rVectorGenerator.EvalSample()[0] + sin(M_PI * c()[1]) * this -> rVectorGenerator.EvalSample()[2]) * exp(-1/8)
+ (cos(M_PI * c()[0]) * this -> rVectorGenerator.EvalSample()[1] + sin(M_PI * c()[0]) * this -> rVectorGenerator.EvalSample()[4])* exp(-1/8));;
}
string Name() const override { return "SparseGridTestPDE";}
};
IStochasticEllipticProblem *
CreateStochasticEllipticProblem(const string &problemName, const Meshes &meshes);
......
#include "TestStochasticCollocation.hpp"
// Test cases for convergence rate
// depending on dimension and level -> Smoothness of problem has to be known
INSTANTIATE_TEST_SUITE_P(
TestStochasticCollocation, TestStochasticCollocationWithoutEpsilon, Values(
TestParams{"SparseGrid2DGeneratorProblem", "FunctionEvaluation",
"DummyPDESolver", 2.513723354063905}
// Todo add more test cases
"DummyPDESolver", 2.513723354063905},
TestParams{"SparseGrid2DPolyDeg2Problem", "FunctionEvaluation",
"DummyPDESolver", 20.0},
TestParams{"SparseGrid2DPD2qptotProblem", "FunctionEvaluation",
"DummyPDESolver", 20.0},
TestParams{"SparseGrid2DHermiteProblem", "FunctionEvaluation",
"DummyPDESolver", 1.0}
));
INSTANTIATE_TEST_SUITE_P(
TestStochasticCollocation, TestStochasticCollocationFailingTests, Values(
TestParams{"SparseGrid2DPD2qptotProblemFail", "FunctionEvaluation",
"DummyPDESolver", 4.2},
TestParams{"SparseGrid2DPD3qptotProblemFail", "FunctionEvaluation",
"DummyPDESolver", 0.9523809523}
));
TEST_P(TestStochasticCollocationWithoutEpsilon, TestSeriellAgainstParallel) {
......@@ -20,6 +36,19 @@ TEST_P(TestStochasticCollocationWithoutEpsilon, TestSeriellAgainstParallel) {
EXPECT_NEAR(scSeriell->aggregate.mean.Q, GetParam().refValue, SC_TEST_TOLERANCE);
}
TEST_P(TestStochasticCollocationFailingTests, Test) {
mout << GetParam() << endl;
mout.StartBlock("Stochastic Collocation seriell");
mout << "Start" << endl;
scSeriell->Method();
mout.EndBlock();
mout << endl;
EXPECT_NE(scSeriell->aggregate.mean.Q, GetParam().refValue);
}
int main(int argc, char **argv) {
return MppTest(
MppTestBuilder(argc, argv).
......
......@@ -67,4 +67,10 @@ public:
TestStochasticCollocation(0.0, 0) {}
};
class TestStochasticCollocationFailingTests : public TestStochasticCollocation {
public:
TestStochasticCollocationFailingTests() :
TestStochasticCollocation(0.0,0) {}
};
#endif //TESTSTOCHASTICCOLLOCATION_HPP
Markdown is supported
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