00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 #include "Sundance.hpp"
00032 #include "SundanceCubicHermite.hpp"
00033 #include "SundanceProblemTesting.hpp"
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047 class SimplePoisson1DTest : public LP1DTestBase
00048 {
00049 public:
00050
00051
00052
00053
00054 SimplePoisson1DTest(const Array<int>& nx)
00055 : LP1DTestBase(nx) {}
00056
00057
00058 std::string name() const {return "SimplePoisson1D";}
00059
00060
00061
00062
00063 Array<int> pExpected() const {return tuple<int>(3);}
00064
00065
00066
00067
00068 Expr exactSoln() const
00069 {
00070 Expr x = coord(0);
00071 const double pi = 4.0*atan(1.0);
00072
00073 return sin(pi*x/2.0);
00074 }
00075
00076
00077
00078
00079
00080 LinearProblem prob(const Mesh& mesh) const
00081 {
00082 const double pi = 4.0*atan(1.0);
00083 CellFilter left = domain().left();
00084
00085 Expr u = new UnknownFunction(new Lagrange(2), "u");
00086 Expr v = new TestFunction(new Lagrange(2), "v");
00087 Expr dx = gradient(1);
00088 Expr x = coord(0);
00089
00090 QuadratureFamily quad = new GaussianQuadrature(4);
00091
00092 Expr eqn = Integral(interior(), (dx*v)*(dx*u), quad)
00093 + Integral(interior(), -0.25*pi*pi*sin(pi*x/2.0)*v, quad);
00094
00095 Expr bc = EssentialBC(left, v*u, quad);
00096
00097 return LinearProblem(mesh, eqn, bc, v, u, vecType());
00098 }
00099
00100
00101
00102
00103 Array<LPTestSpec> specs() const
00104 {
00105
00106
00107 double tol = 0.05;
00108 return tuple(
00109 LPTestSpec("amesos.xml", tol, makeSet<int>(1)),
00110 LPTestSpec("aztec-ifpack.xml", tol),
00111 LPTestSpec("aztec-ml.xml", tol),
00112 LPTestSpec("belos-ml.xml", tol),
00113 LPTestSpec("bicgstab.xml", tol)
00114 );
00115 }
00116 };
00117
00118
00119
00120
00121
00122 class Poisson1DTest : public LP1DTestBase
00123 {
00124 public:
00125
00126 Poisson1DTest(const Array<int>& nx, const BasisFamily& basis)
00127 : LP1DTestBase(nx), basis_(basis) {}
00128
00129
00130 std::string name() const {
00131 TeuchosOStringStream ss;
00132 ss << "Poisson1D(basis=";
00133 basis_.print(ss);
00134 ss << ")";
00135 return ss.str();
00136 }
00137
00138
00139 Array<int> pExpected() const {return tuple<int>(basis_.order()+1);}
00140
00141
00142 Expr exactSoln() const
00143 {
00144 Expr x = coord(0);
00145 const double pi = 4.0*atan(1.0);
00146
00147 return sin(3.0*pi*x/2.0);
00148 }
00149
00150
00151 LinearProblem prob(const Mesh& mesh) const
00152 {
00153 CellFilter left = domain().left();
00154
00155 Expr u = new UnknownFunction(basis_, "u");
00156 Expr v = new TestFunction(basis_, "v");
00157 Expr dx = gradient(1);
00158 Expr x = coord(0);
00159
00160 QuadratureFamily quad = new GaussianQuadrature(2*basis_.order());
00161
00162 const double pi = 4.0*atan(1.0);
00163 Expr eqn = Integral(interior(), (dx*v)*(dx*u), quad)
00164 + Integral(interior(), -2.25*pi*pi*sin(3.0*pi*x/2.0)*v, quad);
00165
00166 Expr bc = EssentialBC(left, v*u, quad);
00167
00168
00169 return LinearProblem(mesh, eqn, bc, v, u, vecType());
00170 }
00171
00172
00173 Array<LPTestSpec> specs() const
00174 {
00175 double tol = 0.05;
00176 return tuple(
00177 LPTestSpec("amesos.xml", tol, makeSet<int>(1)),
00178 LPTestSpec("aztec-ifpack.xml", tol),
00179 LPTestSpec("aztec-ml.xml", tol),
00180 LPTestSpec("belos-ml.xml", tol),
00181 LPTestSpec("bicgstab.xml", tol)
00182 );
00183 }
00184 private:
00185 BasisFamily basis_;
00186 };
00187
00188
00189
00190
00191
00192 class Projection1DTest : public LP1DTestBase
00193 {
00194 public:
00195
00196 Projection1DTest(const Array<int>& nx, const BasisFamily& basis)
00197 : LP1DTestBase(nx), basis_(basis) {}
00198
00199
00200 std::string name() const
00201 {
00202 TeuchosOStringStream ss;
00203 ss << "Projection1D(basis=";
00204 basis_.print(ss);
00205 ss << ")";
00206 return ss.str();
00207 }
00208
00209
00210 Array<int> pExpected() const {return tuple<int>(basis_.order()+1);}
00211
00212
00213 Expr exactSoln() const
00214 {
00215 Expr x = coord(0);
00216 return 10.0*x*exp(-1.0*x);
00217 }
00218
00219
00220 LinearProblem prob(const Mesh& mesh) const
00221 {
00222 Expr u = new UnknownFunction(basis_, "u");
00223 Expr v = new TestFunction(basis_, "v");
00224 Expr x = coord(0);
00225
00226 int p = basis_.order();
00227 QuadratureFamily quad = new GaussianQuadrature(2*p);
00228
00229 Expr ex = exactSoln();
00230 Expr eqn = Integral(interior(), v*(u-ex), quad);
00231 Expr bc;
00232
00233 return LinearProblem(mesh, eqn, bc, v, u, vecType());
00234 }
00235
00236
00237 Array<LPTestSpec> specs() const
00238 {
00239 double tol = 0.05;
00240 return tuple(
00241 LPTestSpec("amesos.xml", tol, makeSet<int>(1)),
00242 LPTestSpec("aztec-ifpack.xml", tol),
00243 LPTestSpec("aztec-ml.xml", tol),
00244 LPTestSpec("bicgstab.xml", tol)
00245 );
00246 }
00247 private:
00248 BasisFamily basis_;
00249 };
00250
00251
00252
00253
00254 class Helmholtz1DTest : public LP1DTestBase
00255 {
00256 public:
00257
00258 Helmholtz1DTest(const Array<int>& nx)
00259 : LP1DTestBase(0.0, atan(1.0), nx) {}
00260
00261
00262 std::string name() const {return "Helmholtz1D";}
00263
00264
00265 Expr exactSoln() const
00266 {
00267 Expr x = coord(0);
00268 return cos(x) + sin(x);
00269 }
00270
00271
00272 Array<int> pExpected() const {return tuple(3);}
00273
00274
00275 LinearProblem prob(const Mesh& mesh) const
00276 {
00277 CellFilter left = domain().left();
00278 CellFilter right = domain().right();
00279
00280 Expr u = new UnknownFunction(new Lagrange(2), "u");
00281 Expr v = new TestFunction(new Lagrange(2), "v");
00282 Expr dx = gradient(1);
00283 Expr x = coord(0);
00284
00285 QuadratureFamily quad = new GaussianQuadrature(4);
00286
00287 Expr eqn = Integral(interior(), (dx*v)*(dx*u) - v*u, quad);
00288 Expr bc = EssentialBC(left, v*(u-cos(x)), quad);
00289
00290 return LinearProblem(mesh, eqn, bc, v, u, vecType());
00291 }
00292
00293
00294 Array<LPTestSpec> specs() const
00295 {
00296 double tol = 0.05;
00297 return tuple(
00298 LPTestSpec("amesos.xml", tol, makeSet<int>(1)),
00299 LPTestSpec("aztec-ifpack.xml", tol),
00300 LPTestSpec("aztec-ml.xml", tol),
00301 LPTestSpec("belos-ml.xml", tol),
00302 LPTestSpec("bicgstab.xml", tol)
00303 );
00304 }
00305 private:
00306 };
00307
00308
00309
00310 class Coupled1DTest : public LP1DTestBase
00311 {
00312 public:
00313
00314 Coupled1DTest(const Array<int>& nx) : LP1DTestBase(nx) {}
00315
00316
00317 std::string name() const {return "Coupled1D";}
00318
00319
00320 Expr exactSoln() const
00321 {
00322 Expr x = coord(0);
00323
00324 Expr x2 = x*x;
00325 Expr x3 = x*x2;
00326
00327 Expr uExact = (1.0/120.0)*x2*x3 - 1.0/36.0 * x3 + 7.0/360.0 * x;
00328 Expr vExact = 1.0/6.0 * x * (x2 - 1.0);
00329
00330 return List(vExact, uExact);
00331 }
00332
00333
00334 Array<int> pExpected() const {return tuple(2, 2);}
00335
00336
00337 LinearProblem prob(const Mesh& mesh) const
00338 {
00339 CellFilter left = domain().left();
00340 CellFilter right = domain().right();
00341
00342
00343 Expr u = new UnknownFunction(new Lagrange(3), "u");
00344 Expr v = new UnknownFunction(new Lagrange(1), "v");
00345 Expr du = new TestFunction(new Lagrange(3), "du");
00346 Expr dv = new TestFunction(new Lagrange(1), "dv");
00347
00348 Expr dx = gradient(1);
00349 Expr x = coord(0);
00350
00351
00352 QuadratureFamily quad = new GaussianQuadrature(10);
00353
00354 Expr eqn = Integral(interior(),
00355 (dx*du)*(dx*u) + du*v + (dx*dv)*(dx*v) + x*dv,
00356 quad);
00357
00358 Expr bc = EssentialBC(left, du*u + dv*v, quad)
00359 + EssentialBC(right, du*u + dv*v, quad);
00360
00361
00362 return LinearProblem(mesh, eqn, bc,
00363 List(dv,du), List(v,u), vecType());
00364
00365 }
00366
00367
00368 Array<LPTestSpec> specs() const
00369 {
00370 double tol = 0.05;
00371 return tuple(
00372 LPTestSpec("amesos.xml", tol, makeSet<int>(1)),
00373 LPTestSpec("aztec-ifpack.xml", tol),
00374 LPTestSpec("belos-ifpack.xml", tol),
00375 LPTestSpec("bicgstab.xml", tol)
00376 );
00377 }
00378 };
00379
00380
00381
00382
00383
00384
00385
00386 int main(int argc, char** argv)
00387 {
00388 try
00389 {
00390 Sundance::init(&argc, &argv);
00391 Tabs::showDepth() = false;
00392 LinearSolveDriver::solveFailureIsFatal() = false;
00393 int np = MPIComm::world().getNProc();
00394
00395 LPTestSuite tests;
00396
00397 Array<int> nx = tuple(8, 16, 24, 32, 40);
00398
00399 tests.registerTest(rcp(new SimplePoisson1DTest(nx)));
00400
00401 for (int p=1; p<=3; p++)
00402 {
00403 if (np==1 || np % 4 == 0)
00404 {
00405 tests.registerTest(rcp(new Poisson1DTest(nx, new Lagrange(p))));
00406 tests.registerTest(rcp(new Poisson1DTest(nx, new Bernstein(p))));
00407 tests.registerTest(rcp(new Projection1DTest(nx, new Lagrange(p))));
00408 tests.registerTest(rcp(new Projection1DTest(nx, new Bernstein(p))));
00409 }
00410 else
00411 {
00412 Out::root() << "skipping 1D tests needing numprocs=1 or 4" << std::endl;
00413 }
00414 }
00415
00416 tests.registerTest(rcp(new Helmholtz1DTest(nx)));
00417
00418
00419 tests.registerTest(rcp(new Coupled1DTest(nx)));
00420
00421 bool pass = tests.run();
00422
00423 Out::root() << "total test status: " << pass << std::endl;
00424
00425 Sundance::passFailTest(pass);
00426 }
00427 catch(std::exception& e)
00428 {
00429 Sundance::handleException(e);
00430 }
00431 Sundance::finalize(); return Sundance::testStatus();
00432 }