51 static double p = 2.0;
54 for (
int i = 0; i < x.
Size(); i++)
56 x(i) = (2*x(i)-xmin(i)-xmax(i))/(xmax(i)-xmin(i));
60 if (x.
Normlp(p) > 0.4 && x.
Normlp(p) < 0.6) {
return 1; }
61 if (x.
Normlp(p) < 0.4 || x.
Normlp(p) > 0.6) {
return 2; }
65 int main(
int argc,
char *argv[])
69 const char *mesh_file =
"../../data/inline-quad.mesh";
74 args.
AddOption(&mesh_file,
"-m",
"--mesh",
75 "Input mesh file to shape materials in.");
76 args.
AddOption(&sd,
"-sd",
"--sub-divisions",
77 "Number of element subdivisions for interface detection.");
78 args.
AddOption(&nclimit,
"-ncl",
"--nc-limit",
79 "Level of hanging nodes allowed (-1 = unlimited).");
80 args.
AddOption(&aniso,
"-a",
"--aniso",
"-i",
"--iso",
81 "Enable anisotropic refinement of quads and hexes.");
87 Mesh mesh(mesh_file, 1, 1);
105 char vishost[] =
"localhost";
108 sol_sock.precision(8);
111 for (
int iter = 0; 1; iter++)
114 for (
int i = 0; i < mesh.
GetNE(); i++)
136 if ((
int)matsum != m*(j+1))
153 int dx = 0, dy = 0, dz = 0;
157 for (
int j = 0; j <= sd; j++)
158 for (
int i = 0; i < sd; i++)
160 dx += abs(mat[j*s + i+1] - mat[j*s + i]);
161 dy += abs(mat[(i+1)*s + j] - mat[i*s + j]);
166 for (
int k = 0; k <= sd; k++)
167 for (
int j = 0; j <= sd; j++)
168 for (
int i = 0; i < sd; i++)
170 dx += abs(mat[(k*s + j)*s + i+1] - mat[(k*s + j)*s + i]);
171 dy += abs(mat[(k*s + i+1)*s + j] - mat[(k*s + i)*s + j]);
172 dz += abs(mat[((i+1)*s + j)*s + k] - mat[(i*s + j)*s + k]);
176 const int tol = mat.Size() / 10;
177 if (dx > tol) { type |= 1; }
178 if (dy > tol) { type |= 2; }
179 if (dz > tol) { type |= 4; }
180 if (!type) { type = 7; }
188 sol_sock <<
"solution\n" << mesh << attr;
189 if (iter == 0 && sdim == 2)
191 sol_sock <<
"keys 'RjlmpppppppppppppA*************'\n";
193 if (iter == 0 && sdim == 3)
195 sol_sock <<
"keys 'YYYYYYYYYXXXXXXXmA********8888888pppttt";
196 if (dim == 3) { sol_sock <<
"iiM"; }
203 cout <<
"Mesh has " << mesh.
GetNE() <<
" elements. \n"
204 <<
"Continue shaping? --> ";
206 if (yn ==
'n' || yn ==
'q') {
break; }
215 for (
int i = 0; i < mesh.
GetNE(); i++)
222 ofstream mesh_ofs(
"shaper.mesh");
223 mesh_ofs.precision(8);
224 mesh.
Print(mesh_ofs);
int GetNPoints() const
Returns the number of the points in the integration rule.
virtual void Print(std::ostream &out=mfem::out) const
Class for an integration rule - an Array of IntegrationPoint.
Class for grid function - Vector with associated FE space.
virtual void Update(bool want_transform=true)
void GetBoundingBox(Vector &min, Vector &max, int ref=2)
RefinedGeometry * Refine(int Geom, int Times, int ETimes=1)
int Size() const
Returns the size of the vector.
int GetNE() const
Returns number of elements.
IntegrationPoint & IntPoint(int i)
Returns a reference to the i-th integration point.
int material(Vector &x, Vector &xmin, Vector &xmax)
double Normlp(double p) const
Returns the l_p norm of the vector.
int Append(const T &el)
Append element to array, resize if necessary.
void SetCurvature(int order, bool discont=false, int space_dim=-1, int ordering=1)
int MeshGenerator()
Get the mesh generator/type.
GeometryRefiner GlobGeometryRefiner
void PrintUsage(std::ostream &out) const
int SpaceDimension() const
void AddOption(bool *var, const char *enable_short_name, const char *enable_long_name, const char *disable_short_name, const char *disable_long_name, const char *description, bool required=false)
int GetElementBaseGeometry(int i=0) const
NURBSExtension * NURBSext
Optional NURBS mesh extension.
void GetElementTransformation(int i, IsoparametricTransformation *ElTr)
void PrintOptions(std::ostream &out) const
void SetAttribute(int i, int attr)
Set the attribute of element i.
void GeneralRefinement(const Array< Refinement > &refinements, int nonconforming=-1, int nc_limit=0)
Arbitrary order "L2-conforming" discontinuous finite elements.