Commit ce676806 authored by Kévin Moreau's avatar Kévin Moreau

update mode atomic test case

parent 0d20941f
#include "FastMarchingInterface.h"
#include "TLSAlgorithm.h"
#include "xAlgorithm.h"
#include "xLevelSet.h"
#include "xMesh.h"
......@@ -25,19 +23,25 @@ namespace lalg {
using namespace xfem;
using namespace lalg;
using namespace xtls;
using Trellis_Util::mPoint;
using AOMD::mEntity;
using AOMD::mVertex;
int main(int argc, char* argv[]) {
if(argc!=3) {
std::cout<<"First argument: .msh file; Second argument: a boolean to say if bnd_in is to do."<<std::endl;
std::abort();
}
xMesh mesh(argv[1]);
// domain creation
xPlane p1(mPoint(0.1,0.1,0.), xVector(1.,0.,0.));
xPlane p2(mPoint(0.205,0.205,0.), xVector(-1.,0.,0.));
xLevelSet level_set(&mesh, xUnion(p1, p2));
xGrownHalfLine grown(mPoint(0.5, 0.5, 0.), xVector(1.,0.,0.), 0.26);
xLevelSet level_set(&mesh, xCompl(grown));
// xPlane p1(mPoint(0.1,0.1,0.), xVector(1.,0.,0.));
// xPlane p2(mPoint(0.205,0.205,0.), xVector(-1.,0.,0.));
// xLevelSet level_set(&mesh, xUnion(p1, p2));
xPhysSurfByTagging phys_surf(level_set);
xEntityFilter filter=std::bind1st(std::mem_fun(&xPhysSurfByTagging::strictOut), &phys_surf);
xSubMesh domain("domain", mesh);
......@@ -69,21 +73,22 @@ int main(int argc, char* argv[]) {
}
}
xRegion all(&mesh);
xRegion region(&domain);
std::cout<<"region "<<region.size(2)<<std::endl;
std::cout<<"bnd "<<bnd.size(0)<<std::endl;
std::cout<<"interior "<<interior.size(0)<<std::endl;
level_set.load(p1);
// level_set.load(p1);
// phi field
xDoubleManager double_manager;
xSpacePolynomialLagrange space("space", xSpace::SCALAR, 1);
xField field(&double_manager, space);
DeclareInterpolation(field, xValueCreator<xValueDouble>(), region.begin(), region.end());
DeclareState(field, xStateDofCreator<>(double_manager, "dofs"), region.begin(), region.end());
for(xIter it=region.begin(0); it!=region.end(0); ++it) {
DeclareInterpolation(field, xValueCreator<xValueDouble>(), all.begin(), all.end());
DeclareState(field, xStateDofCreator<>(double_manager, "dofs"), all.begin(), all.end());
for(xIter it=all.begin(0); it!=all.end(0); ++it) {
xFiniteElementKeysOnly keys;
keys.setKeys(*it, field.begin(), field.end());
std::vector<xValue<double>*> vals;
......@@ -95,7 +100,7 @@ int main(int argc, char* argv[]) {
xEvalField<xIdentity<double> > eval_field(field);
xIntegrationRuleBasic integ_rule(1);
xExportGmshAscii pexport;
Export(eval_field, pexport, "phi_given", integ_rule, region.begin(), region.end());
Export(eval_field, pexport, "phi_given", integ_rule, all.begin(), all.end());
std::function<double (const AOMD::mVertex&)> f_func=[](const AOMD::mVertex&){ return 1.; };
......@@ -106,13 +111,20 @@ int main(int argc, char* argv[]) {
region,
bnd.begin(0), bnd.end(0),
interior.begin(0), interior.end(0),
std::back_inserter(bnd_in_vec));
std::back_inserter(bnd_in_vec), 1.e-6);
Export(eval_field, pexport, "phi", integ_rule, region.begin(), region.end());
Export(eval_field, pexport, "phi", integ_rule, all.begin(), all.end());
xSubMesh bnd_in("bnd_in", mesh);
for(auto v: bnd_in_vec) {
bnd_in.add(const_cast<AOMD::mEntity*>(v));
if(atoi(argv[2])>0) {
for(auto v: bnd_in_vec) {
bnd_in.add(const_cast<AOMD::mEntity*>(v));
}
}
else {
for(auto it=bnd.begin(0); it!=bnd.end(0); ++it) {
bnd_in.add(*it);
}
}
xSubMesh interior_in("interior_in", mesh);
xAcceptInSubMesh is_in_bnd_in(bnd_in);
......@@ -133,12 +145,20 @@ int main(int argc, char* argv[]) {
region,
bnd_in.begin(0), bnd_in.end(0),
interior_in.begin(0), interior_in.end(0),
coeffs);
coeffs, 1.e-6);
std::cout<<"nb modes "<<bnd_in.size(0)<<std::endl;
// modes export
for(int i=0; i<bnd_in.size(0); ++i) {
for(xIter it=all.begin(0); it!=all.end(0); ++it) {
xFiniteElementKeysOnly keys;
keys.setKeys(*it, field.begin(), field.end());
std::vector<xValue<double>*> vals;
vals.reserve(keys.sizeKey());
double_manager.getValPtr(keys.beginKey(), keys.endKey(), vals);
vals[0]->setVal(0.);
}
int j=0;
for(xIter it=region.begin(0); it!=region.end(0); ++it, ++j) {
xFiniteElementKeysOnly keys;
......@@ -150,7 +170,7 @@ int main(int argc, char* argv[]) {
}
std::ostringstream oss;
oss<<i;
Export(eval_field, pexport, "mode_"+oss.str(), integ_rule, region.begin(), region.end());
Export(eval_field, pexport, "mode_"+oss.str(), integ_rule, all.begin(), all.end());
}
int j=0;
......@@ -167,7 +187,7 @@ int main(int argc, char* argv[]) {
}
vals[0]->setVal(sum);
}
Export(eval_field, pexport, "sum", integ_rule, region.begin(), region.end());
Export(eval_field, pexport, "sum", integ_rule, all.begin(), all.end());
return 0;
}
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