# To specify which files are to be deleted by "make clean" (apart from
# the usual ones: object files, executables, backups, etc), fill in the
# following list
-delete-files = gnuplot* *.eps
+delete-files = *gnuplot *.eps *ucd *history
#include <basic/function.h>
#include <basic/data_io.h>
#include <fe/fe_lib.h>
+#include <fe/fe_lib.criss_cross.h>
#include <fe/quadrature_lib.h>
#include <numerics/base.h>
#include <numerics/assembler.h>
create_new ();
cout << "Refinement level = " << level
- << ", using elements of order " << order
- << endl;
+ << ", using elements of type <";
+ switch (order)
+ {
+ case 0:
+ cout << "criss-cross";
+ break;
+ default:
+ cout << "Lagrange-" << order;
+ break;
+ };
+ cout << ">" << endl;
cout << " Making grid... ";
tria->create_hyper_ball ();
Quadrature<dim> *quadrature;
Quadrature<dim-1> *boundary_quadrature;
switch (order) {
+ case 0:
+ fe = new FECrissCross<dim>();
+ quadrature = new QCrissCross1<dim>();
+ boundary_quadrature = new QGauss2<dim-1>();
+ break;
case 1:
fe = new FELinear<dim>();
quadrature = new QGauss3<dim>();
h1_seminorm_error_per_dof);
dof->distribute_cell_to_dof_vector (h1_error_per_cell, h1_error_per_dof);
- dVector projected_solution;
- ConstraintMatrix constraints;
- constraints.close ();
- VectorTools<dim>::project (*dof, constraints, *fe,
- StraightBoundary<dim>(), *quadrature,
- sol, projected_solution, false,
- *boundary_quadrature);
- cout << " Calculating L2 error of projected solution... ";
- VectorTools<dim>::integrate_difference (*dof_handler,
- projected_solution, sol,
- l2_error_per_cell,
- *quadrature, *fe, L2_norm);
- cout << l2_error_per_cell.l2_norm() << endl;
+// dVector projected_solution;
+// ConstraintMatrix constraints;
+// constraints.close ();
+// VectorTools<dim>::project (*dof, constraints, *fe,
+// StraightBoundary<dim>(), *quadrature,
+// sol, projected_solution, false,
+// *boundary_quadrature);
+// cout << " Calculating L2 error of projected solution... ";
+// VectorTools<dim>::integrate_difference (*dof_handler,
+// projected_solution, sol,
+// l2_error_per_cell,
+// *quadrature, *fe, L2_norm);
+// cout << l2_error_per_cell.l2_norm() << endl;
string filename;
DataOut<dim> out;
ofstream o(filename.c_str());
fill_data (out);
- out.add_data_vector (projected_solution, "projected u");
out.add_data_vector (l1_error_per_dof, "L1-Error");
out.add_data_vector (l2_error_per_dof, "L2-Error");
out.add_data_vector (linfty_error_per_dof, "Linfty-Error");
int main () {
- for (unsigned int order=1; order<5; ++order)
+ for (unsigned int order=0; order<5; ++order)
{
PoissonProblem<2> problem (order);
string filename;
switch (order)
{
+ case 0:
+ filename = "criss_cross";
+ break;
case 1:
filename = "linear";
break;
+set term postscript eps
set xlabel "Number of degrees of freedom"
-set ylabel "Error"
set data style linespoints
set logscale xy
-set term postscript eps
+
+
+set ylabel "Error"
+
+set output "criss-cross.eps"
+
+plot "criss-cross.history" using 1:2 title "L1 error","criss-cross.history" using 1:3 title "L2 error","criss-cross.history" using 1:4 title "Linfty error","criss-cross.history" using 1:5 title "H1 seminorm error","criss-cross.history" using 1:6 title "H1 error"
+
+
+
set output "linear.eps"
plot "linear.history" using 1:2 title "L1 error","linear.history" using 1:3 title "L2 error","linear.history" using 1:4 title "Linfty error","linear.history" using 1:5 title "H1 seminorm error","linear.history" using 1:6 title "H1 error"
-set term postscript eps
set output "quadratic.eps"
plot "quadratic.history" using 1:2 title "L1 error","quadratic.history" using 1:3 title "L2 error","quadratic.history" using 1:4 title "Linfty error","quadratic.history" using 1:5 title "H1 seminorm error","quadratic.history" using 1:6 title "H1 error"
-set term postscript eps
set output "cubic.eps"
plot "cubic.history" using 1:2 title "L1 error","cubic.history" using 1:3 title "L2 error","cubic.history" using 1:4 title "Linfty error","cubic.history" using 1:5 title "H1 seminorm error","cubic.history" using 1:6 title "H1 error"
-set term postscript eps
set output "quartic.eps"
plot "quartic.history" using 1:2 title "L1 error","quartic.history" using 1:3 title "L2 error","quartic.history" using 1:4 title "Linfty error","quartic.history" using 1:5 title "H1 seminorm error","quartic.history" using 1:6 title "H1 error"
-set term postscript eps
set output "l2error.eps"
set ylabel "L2-error"
-plot "linear.history" using 1:3 title "Linear elements", "quadratic.history" using 1:3 title "Quadratic elements", "cubic.history" using 1:3 title "Cubic elements", "quartic.history" using 1:3 title "Quartic elements"
+plot "criss-cross.history" using 1:3 title "Criss-cross elements", "linear.history" using 1:3 title "Linear elements", "quadratic.history" using 1:3 title "Quadratic elements", "cubic.history" using 1:3 title "Cubic elements", "quartic.history" using 1:3 title "Quartic elements"
-set term postscript eps
set output "h1error.eps"
set ylabel "H1-error"
-plot "linear.history" using 1:6 title "Linear elements", "quadratic.history" using 1:6 title "Quadratic elements", "cubic.history" using 1:6 title "Cubic elements", "quartic.history" using 1:6 title "Quartic elements"
+plot "criss-cross.history" using 1:6 title "Criss-cross elements", "linear.history" using 1:6 title "Linear elements", "quadratic.history" using 1:6 title "Quadratic elements", "cubic.history" using 1:6 title "Cubic elements", "quartic.history" using 1:6 title "Quartic elements"
# To specify which files are to be deleted by "make clean" (apart from
# the usual ones: object files, executables, backups, etc), fill in the
# following list
-delete-files = gnuplot* *.eps
+delete-files = *gnuplot *.eps *ucd *history
#include <basic/function.h>
#include <basic/data_io.h>
#include <fe/fe_lib.h>
+#include <fe/fe_lib.criss_cross.h>
#include <fe/quadrature_lib.h>
#include <numerics/base.h>
#include <numerics/assembler.h>
create_new ();
cout << "Refinement level = " << level
- << ", using elements of order " << order
- << endl;
+ << ", using elements of type <";
+ switch (order)
+ {
+ case 0:
+ cout << "criss-cross";
+ break;
+ default:
+ cout << "Lagrange-" << order;
+ break;
+ };
+ cout << ">" << endl;
cout << " Making grid... ";
tria->create_hyper_ball ();
Quadrature<dim> *quadrature;
Quadrature<dim-1> *boundary_quadrature;
switch (order) {
+ case 0:
+ fe = new FECrissCross<dim>();
+ quadrature = new QCrissCross1<dim>();
+ boundary_quadrature = new QGauss2<dim-1>();
+ break;
case 1:
fe = new FELinear<dim>();
quadrature = new QGauss3<dim>();
h1_seminorm_error_per_dof);
dof->distribute_cell_to_dof_vector (h1_error_per_cell, h1_error_per_dof);
- dVector projected_solution;
- ConstraintMatrix constraints;
- constraints.close ();
- VectorTools<dim>::project (*dof, constraints, *fe,
- StraightBoundary<dim>(), *quadrature,
- sol, projected_solution, false,
- *boundary_quadrature);
- cout << " Calculating L2 error of projected solution... ";
- VectorTools<dim>::integrate_difference (*dof_handler,
- projected_solution, sol,
- l2_error_per_cell,
- *quadrature, *fe, L2_norm);
- cout << l2_error_per_cell.l2_norm() << endl;
+// dVector projected_solution;
+// ConstraintMatrix constraints;
+// constraints.close ();
+// VectorTools<dim>::project (*dof, constraints, *fe,
+// StraightBoundary<dim>(), *quadrature,
+// sol, projected_solution, false,
+// *boundary_quadrature);
+// cout << " Calculating L2 error of projected solution... ";
+// VectorTools<dim>::integrate_difference (*dof_handler,
+// projected_solution, sol,
+// l2_error_per_cell,
+// *quadrature, *fe, L2_norm);
+// cout << l2_error_per_cell.l2_norm() << endl;
string filename;
DataOut<dim> out;
ofstream o(filename.c_str());
fill_data (out);
- out.add_data_vector (projected_solution, "projected u");
out.add_data_vector (l1_error_per_dof, "L1-Error");
out.add_data_vector (l2_error_per_dof, "L2-Error");
out.add_data_vector (linfty_error_per_dof, "Linfty-Error");
int main () {
- for (unsigned int order=1; order<5; ++order)
+ for (unsigned int order=0; order<5; ++order)
{
PoissonProblem<2> problem (order);
string filename;
switch (order)
{
+ case 0:
+ filename = "criss_cross";
+ break;
case 1:
filename = "linear";
break;
+set term postscript eps
set xlabel "Number of degrees of freedom"
-set ylabel "Error"
set data style linespoints
set logscale xy
-set term postscript eps
+
+
+set ylabel "Error"
+
+set output "criss-cross.eps"
+
+plot "criss-cross.history" using 1:2 title "L1 error","criss-cross.history" using 1:3 title "L2 error","criss-cross.history" using 1:4 title "Linfty error","criss-cross.history" using 1:5 title "H1 seminorm error","criss-cross.history" using 1:6 title "H1 error"
+
+
+
set output "linear.eps"
plot "linear.history" using 1:2 title "L1 error","linear.history" using 1:3 title "L2 error","linear.history" using 1:4 title "Linfty error","linear.history" using 1:5 title "H1 seminorm error","linear.history" using 1:6 title "H1 error"
-set term postscript eps
set output "quadratic.eps"
plot "quadratic.history" using 1:2 title "L1 error","quadratic.history" using 1:3 title "L2 error","quadratic.history" using 1:4 title "Linfty error","quadratic.history" using 1:5 title "H1 seminorm error","quadratic.history" using 1:6 title "H1 error"
-set term postscript eps
set output "cubic.eps"
plot "cubic.history" using 1:2 title "L1 error","cubic.history" using 1:3 title "L2 error","cubic.history" using 1:4 title "Linfty error","cubic.history" using 1:5 title "H1 seminorm error","cubic.history" using 1:6 title "H1 error"
-set term postscript eps
set output "quartic.eps"
plot "quartic.history" using 1:2 title "L1 error","quartic.history" using 1:3 title "L2 error","quartic.history" using 1:4 title "Linfty error","quartic.history" using 1:5 title "H1 seminorm error","quartic.history" using 1:6 title "H1 error"
-set term postscript eps
set output "l2error.eps"
set ylabel "L2-error"
-plot "linear.history" using 1:3 title "Linear elements", "quadratic.history" using 1:3 title "Quadratic elements", "cubic.history" using 1:3 title "Cubic elements", "quartic.history" using 1:3 title "Quartic elements"
+plot "criss-cross.history" using 1:3 title "Criss-cross elements", "linear.history" using 1:3 title "Linear elements", "quadratic.history" using 1:3 title "Quadratic elements", "cubic.history" using 1:3 title "Cubic elements", "quartic.history" using 1:3 title "Quartic elements"
-set term postscript eps
set output "h1error.eps"
set ylabel "H1-error"
-plot "linear.history" using 1:6 title "Linear elements", "quadratic.history" using 1:6 title "Quadratic elements", "cubic.history" using 1:6 title "Cubic elements", "quartic.history" using 1:6 title "Quartic elements"
+plot "criss-cross.history" using 1:6 title "Criss-cross elements", "linear.history" using 1:6 title "Linear elements", "quadratic.history" using 1:6 title "Quadratic elements", "cubic.history" using 1:6 title "Cubic elements", "quartic.history" using 1:6 title "Quartic elements"