Skip to content
Snippets Groups Projects
Commit 8b26c10d authored by Boris Martin's avatar Boris Martin
Browse files

cleanup, documentation

parent bf173fda
No related branches found
No related tags found
1 merge request!17New tutorial case: overlapping DDM (Dirichlet) for Laplace equation
......@@ -29,22 +29,27 @@ int main(int argc, char **argv)
GmshDdm gmshDdm(argc, argv);
unsigned int nDom = 2;
gmshDdm.userDefinedParameter(nDom, "nDom");
double lc = 0.05;
int order = 1;
double lc = 0.05;
std::string gauss = "Gauss10";
gmshDdm.userDefinedParameter(lc, "lc");
gmshDdm.userDefinedParameter(nDom, "nDom");
gmshDdm.userDefinedParameter(order, "order");
std::string gauss = "Gauss10";
gmshDdm.userDefinedParameter(lc, "lc");
gmshDdm.userDefinedParameter(gauss, "gauss");
gmsh::open("../domain.geo");
gmsh::open("../domain.geo"); // Assume we run from the "build" subfolder
gmsh::model::mesh::generate();
Domain left(2, 11), center(2, 12), right(2, 13);
Domain gammaLeft(1, 14), interface1(1, 16), interface2(1, 17), gammaRight(1, 15);
// Domains are partitions of the x-axis.
Domain left("left"); // 0 to 1.2
Domain center("center"); // 0.8 to 1.2
Domain right("right"); // 0.8 to 2
Domain gammaLeft("gammaLeft"); // 0
Domain interface1("interface1"); // 0.8
Domain interface2("interface2"); // 1.2
Domain gammaRight("gammaRight"); // 2.0
Subdomain omega(nDom);
Interface sigma(nDom);
......@@ -59,78 +64,54 @@ int main(int argc, char **argv)
topology[1] = {0};
// Create DDM formulation
// And regular BCs
// And BCs on the subdomain problems
gmshddm::problem::Formulation< std::complex< double > > formulation("Laplace", topology);
SubdomainField< std::complex< double >, Form::Form0 > u("u", omega, FunctionSpaceTypeForm0::HierarchicalH1, order);
// Toy problem with solution u(x,y) = x
u(0).addConstraint(gammaLeft, 0);
u(1).addConstraint(gammaRight, 1);
u(1).addConstraint(gammaRight, 2);
auto trueSolution = x<std::complex< double >>();
InterfaceField< std::complex< double >, Form::Form0 > g("g", sigma, FunctionSpaceTypeForm0::HierarchicalH1, order);
InterfaceField< std::complex< double >, Form::Form0 > g("g", sigma, FunctionSpaceTypeForm0::HierarchicalH1, order);
formulation.addInterfaceField(g);
//u(0).addConstraint(sigma(0,1), g(1,0));
// First domain
formulation(0).integral(grad(dof(u(0))), grad(tf(u(0))), omega(0), gauss);
//formulation(0).integral((-g(1, 0)), tf(u(0)), sigma(0, 1), gauss);
formulation(0).integral((-g(1, 0)), tf(u(0)), sigma(0, 1), gauss);
formulation(0).integral(dof(u(0)), tf(u(0)), sigma(0, 1), gauss);
//formulation(0).integral(formulation.artificialSource(-g(1, 0)), tf(u(0)), sigma(0, 1), gauss);
formulation(0).integral((-g(1, 0)), tf(u(0)), interface2, gauss);
formulation(0).integral(dof(u(0)), tf(u(0)), interface2, gauss);
// First interface
// First interface (What 0 sends to 1)
formulation(0, 1).integral(-dof(g(0, 1)), tf(g(0, 1)), sigma(0, 1), gauss);
formulation(0, 1).integral(-(g(1, 0)), tf(g(0, 1)), sigma(0, 1), gauss);
formulation(0, 1).integral((2*u(0)), tf(g(0, 1)), sigma(0, 1), gauss);
//formulation(0, 1).integral(dof(u(0)), tf(g(0, 1)), sigma(0, 1), gauss);
//formulation(0, 1).integral(formulation.artificialSource(g(1, 0)), tf(g(0, 1)), sigma(0, 1), gauss);
// Second domain
formulation(1).integral(grad(dof(u(1))), grad(tf(u(1))), omega(1), gauss);
//formulation(1).integral(formulation.artificialSource(-g(0, 1)), tf(u(1)), sigma(1, 0), gauss);
formulation(1).integral((-g(0, 1)), tf(u(1)), sigma(1, 0), gauss);
formulation(1).integral(dof(u(1)), tf(u(1)), sigma(1, 0), gauss);
formulation(1).integral((-g(0, 1)), tf(u(1)), interface1, gauss);
formulation(1).integral(dof(u(1)), tf(u(1)), interface1, gauss);
// Second interface
formulation(1, 0).integral(-dof(g(1, 0)), tf(g(1, 0)), sigma(1, 0), gauss);
formulation(1, 0).integral(-(g(0, 1)), tf(g(1, 0)), sigma(1, 0), gauss);
formulation(1, 0).integral((2*u(1)), tf(g(1, 0)), sigma(1, 0), gauss);
//formulation(1, 0).integral(formulation.artificialSource(g(0, 1)), tf(g(1, 0)), sigma(1, 0), gauss);
std::cout << "End formulation g\n";
formulation.pre();
int iterMax = 1000;
gmshDdm.userDefinedParameter(iterMax, "iter");
formulation.solve("gmres", 1e-6, iterMax);
for(unsigned int i = 0; i < nDom; ++i) {
if(isItMySubdomain(i)) {
save(u(i), omega(i), "u_" + std::to_string(i));
save(u(i) - trueSolution, omega(i), "err_" + std::to_string(i));
}
}
save(g(0,1));
save(g(1,0));
/*
Formulation< std::complex< double > > formulation("demoFloatingPotential");
Field<std::complex<double>, Form::Form0> v("v", left|center|interface2|interface1, FunctionSpaceTypeForm0::Lagrange);
v.addConstraint(interface2, 0);
formulation.integral(grad(dof(v)), grad(tf(v)), left|center|interface2|interface1, gauss);
formulation.integral(-1.0, tf(v), gammaLeft, gauss);
formulation.pre();
formulation.assemble();
formulation.solve();
save(v);*/
return 0;
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment