From b6b80eb955bf63445cbae2b4a7ae78676680c721 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Fri, 28 Aug 1998 09:08:04 +0000 Subject: [PATCH] Add a testcase for the multigrid components. git-svn-id: https://svn.dealii.org/trunk@531 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/Attic/examples/Makefile | 2 +- .../deal.II/Attic/examples/multigrid/Makefile | 118 ++++ .../deal.II/Attic/examples/multigrid/make_ps | 52 ++ tests/big-tests/Makefile | 2 +- tests/big-tests/multigrid/Makefile | 118 ++++ tests/big-tests/multigrid/convergence.cc | 511 ++++++++++++++++++ tests/big-tests/multigrid/make_ps | 52 ++ 7 files changed, 853 insertions(+), 2 deletions(-) create mode 100644 deal.II/deal.II/Attic/examples/multigrid/Makefile create mode 100644 deal.II/deal.II/Attic/examples/multigrid/make_ps create mode 100644 tests/big-tests/multigrid/Makefile create mode 100644 tests/big-tests/multigrid/convergence.cc create mode 100644 tests/big-tests/multigrid/make_ps diff --git a/deal.II/deal.II/Attic/examples/Makefile b/deal.II/deal.II/Attic/examples/Makefile index 3c2149f907..1a75c07be4 100644 --- a/deal.II/deal.II/Attic/examples/Makefile +++ b/deal.II/deal.II/Attic/examples/Makefile @@ -3,7 +3,7 @@ # list the directories we want to visit -subdirs = grid/ dof/ poisson/ convergence/ error-estimation/ +subdirs = grid/ dof/ poisson/ convergence/ error-estimation/ multigrid/ # define lists of targets: for each directory we produce a target name # for compilation, running and cleaning by appending the action to diff --git a/deal.II/deal.II/Attic/examples/multigrid/Makefile b/deal.II/deal.II/Attic/examples/multigrid/Makefile new file mode 100644 index 0000000000..022319dc4e --- /dev/null +++ b/deal.II/deal.II/Attic/examples/multigrid/Makefile @@ -0,0 +1,118 @@ +# $Id$ +# Copyright W. Bangerth, University of Heidelberg, 1998 + +# Template for makefiles for the examples subdirectory. In principle, +# everything should be done automatically if you set the target file +# here correctly: +target = convergence + +# All dependencies between files should be updated by the included +# file Makefile.dep if necessary. Object files are compiled into +# the archives ./Obj.a and ./Obj.g.a. By default, the debug version +# is used to link. It you don't like that, change the following +# variable to "off" +debug-mode = off + +# If you want your program to be linked with extra object or library +# files, specify them here: +user-libs = + +# To run the program, use "make run"; to give parameters to the program, +# give the parameters to the following variable: +run-parameters = + +# To execute additional action apart from running the program, fill +# in this list: +additional-run-action = gnuplot make_ps + +# 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 *ucd *history + + + + +############################################################################### +# Internals + +include ../../Make.global_options + + +# get lists of files we need +cc-files = $(wildcard *.cc) +o-files = $(cc-files:.cc=.o) +go-files = $(cc-files:.cc=.go) +h-files = $(wildcard *.h) +lib-h-files = $(wildcard ../../include/*/*.h) + +# list of libraries needed to link with +libs = ./Obj.a ../../lib/libdeal_II.a ../../../lac/lib/liblac.a +libs.g = ./Obj.g.a ../../lib/libdeal_II.g.a ../../../lac/lib/liblac.g.a + + +# check whether we use debug mode or not +ifeq ($(debug-mode),on) +libraries = $(libs.g) +flags = $(CXXFLAGS.g) +endif + +ifeq ($(debug-mode),off) +libraries = $(libs) +flags = $(CXXFLAGS) +endif + + + +# make rule for the target +$(target) : $(libraries) + @echo ============================ Linking $@ + @$(CXX) $(flags) -o $@ $(libraries) $(user-libs) + +# rule how to run the program +run: $(target) + $(target) $(run-parameters) + $(additional-run-action) + + +# rule to make object files +%.go : %.cc + @echo ============================ Compiling with debugging information: $< + @echo $(CXX) ... -c $< -o $@ + @$(CXX) $(CXXFLAGS.g) -c $< -o $@ +%.o : %.cc + @echo ============================ Compiling with optimization: $< + @echo $(CXX) ... -c $< -o $@ + @$(CXX) $(CXXFLAGS) -c $< -o $@ + + +# rules which files the libraries depend upon +Obj.a: ./Obj.a($(o-files)) +Obj.g.a: ./Obj.g.a($(go-files)) + + +clean: + rm -f *.o *.go *~ Makefile.dep Obj.a Obj.g.a $(target) $(delete-files) + + + +.PHONY: clean + + +#Rule to generate the dependency file. This file is +#automagically remade whenever needed, i.e. whenever +#one of the cc-/h-files changed. Make detects whether +#to remake this file upon inclusion at the bottom +#of this file. +# +#use perl to generate rules for the .go files as well +#as to make rules not for tria.o and the like, but +#rather for libnumerics.a(tria.o) +Makefile.dep: $(cc-files) $(h-files) $(lib-h-files) + @echo ============================ Remaking Makefile + @perl ../../Make_dep.pl ./Obj $(INCLUDE) $(cc-files) \ + > Makefile.dep + + +include Makefile.dep + diff --git a/deal.II/deal.II/Attic/examples/multigrid/make_ps b/deal.II/deal.II/Attic/examples/multigrid/make_ps new file mode 100644 index 0000000000..c2df845aec --- /dev/null +++ b/deal.II/deal.II/Attic/examples/multigrid/make_ps @@ -0,0 +1,52 @@ +set term postscript eps +set xlabel "Number of degrees of freedom" +set data style linespoints +set logscale xy + + + +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 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 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 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 output "l2error.eps" +set ylabel "L2-error" + +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 output "h1error.eps" +set ylabel "H1-error" + +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" + + diff --git a/tests/big-tests/Makefile b/tests/big-tests/Makefile index 3c2149f907..1a75c07be4 100644 --- a/tests/big-tests/Makefile +++ b/tests/big-tests/Makefile @@ -3,7 +3,7 @@ # list the directories we want to visit -subdirs = grid/ dof/ poisson/ convergence/ error-estimation/ +subdirs = grid/ dof/ poisson/ convergence/ error-estimation/ multigrid/ # define lists of targets: for each directory we produce a target name # for compilation, running and cleaning by appending the action to diff --git a/tests/big-tests/multigrid/Makefile b/tests/big-tests/multigrid/Makefile new file mode 100644 index 0000000000..022319dc4e --- /dev/null +++ b/tests/big-tests/multigrid/Makefile @@ -0,0 +1,118 @@ +# $Id$ +# Copyright W. Bangerth, University of Heidelberg, 1998 + +# Template for makefiles for the examples subdirectory. In principle, +# everything should be done automatically if you set the target file +# here correctly: +target = convergence + +# All dependencies between files should be updated by the included +# file Makefile.dep if necessary. Object files are compiled into +# the archives ./Obj.a and ./Obj.g.a. By default, the debug version +# is used to link. It you don't like that, change the following +# variable to "off" +debug-mode = off + +# If you want your program to be linked with extra object or library +# files, specify them here: +user-libs = + +# To run the program, use "make run"; to give parameters to the program, +# give the parameters to the following variable: +run-parameters = + +# To execute additional action apart from running the program, fill +# in this list: +additional-run-action = gnuplot make_ps + +# 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 *ucd *history + + + + +############################################################################### +# Internals + +include ../../Make.global_options + + +# get lists of files we need +cc-files = $(wildcard *.cc) +o-files = $(cc-files:.cc=.o) +go-files = $(cc-files:.cc=.go) +h-files = $(wildcard *.h) +lib-h-files = $(wildcard ../../include/*/*.h) + +# list of libraries needed to link with +libs = ./Obj.a ../../lib/libdeal_II.a ../../../lac/lib/liblac.a +libs.g = ./Obj.g.a ../../lib/libdeal_II.g.a ../../../lac/lib/liblac.g.a + + +# check whether we use debug mode or not +ifeq ($(debug-mode),on) +libraries = $(libs.g) +flags = $(CXXFLAGS.g) +endif + +ifeq ($(debug-mode),off) +libraries = $(libs) +flags = $(CXXFLAGS) +endif + + + +# make rule for the target +$(target) : $(libraries) + @echo ============================ Linking $@ + @$(CXX) $(flags) -o $@ $(libraries) $(user-libs) + +# rule how to run the program +run: $(target) + $(target) $(run-parameters) + $(additional-run-action) + + +# rule to make object files +%.go : %.cc + @echo ============================ Compiling with debugging information: $< + @echo $(CXX) ... -c $< -o $@ + @$(CXX) $(CXXFLAGS.g) -c $< -o $@ +%.o : %.cc + @echo ============================ Compiling with optimization: $< + @echo $(CXX) ... -c $< -o $@ + @$(CXX) $(CXXFLAGS) -c $< -o $@ + + +# rules which files the libraries depend upon +Obj.a: ./Obj.a($(o-files)) +Obj.g.a: ./Obj.g.a($(go-files)) + + +clean: + rm -f *.o *.go *~ Makefile.dep Obj.a Obj.g.a $(target) $(delete-files) + + + +.PHONY: clean + + +#Rule to generate the dependency file. This file is +#automagically remade whenever needed, i.e. whenever +#one of the cc-/h-files changed. Make detects whether +#to remake this file upon inclusion at the bottom +#of this file. +# +#use perl to generate rules for the .go files as well +#as to make rules not for tria.o and the like, but +#rather for libnumerics.a(tria.o) +Makefile.dep: $(cc-files) $(h-files) $(lib-h-files) + @echo ============================ Remaking Makefile + @perl ../../Make_dep.pl ./Obj $(INCLUDE) $(cc-files) \ + > Makefile.dep + + +include Makefile.dep + diff --git a/tests/big-tests/multigrid/convergence.cc b/tests/big-tests/multigrid/convergence.cc new file mode 100644 index 0000000000..0aa8f561f7 --- /dev/null +++ b/tests/big-tests/multigrid/convergence.cc @@ -0,0 +1,511 @@ +/* $Id$ */ +/* Copyright W. Bangerth, University of Heidelberg, 1998 */ + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include + + + + + +template +class PoissonEquation : public Equation { + public: + PoissonEquation (const Function &rhs) : + Equation(1), + right_hand_side (rhs) {}; + + virtual void assemble (dFMatrix &cell_matrix, + dVector &rhs, + const FEValues &fe_values, + const Triangulation::cell_iterator &cell) const; + virtual void assemble (dFMatrix &cell_matrix, + const FEValues &fe_values, + const Triangulation::cell_iterator &cell) const; + virtual void assemble (dVector &rhs, + const FEValues &fe_values, + const Triangulation::cell_iterator &cell) const; + protected: + const Function &right_hand_side; +}; + + + + + + +template +class PoissonProblem : public ProblemBase { + public: + PoissonProblem (unsigned int order); + + void clear (); + void create_new (); + int run (unsigned int level); + void print_history (string filename) const; + + protected: + Triangulation *tria; + DoFHandler *dof; + + Function *rhs; + Function *boundary_values; + + vector l1_error, l2_error, linfty_error, h1_seminorm_error, h1_error; + vector n_dofs; + + unsigned int order; +}; + + + + + +/** + Right hand side constructed such that the exact solution is + $sin(2 pi x) + sin(2 pi y)$ + */ +template +class RHSPoly : public Function { + public: + /** + * Return the value of the function + * at the given point. + */ + virtual double operator () (const Point &p) const; +}; + + + +template +class Solution : public Function { + public: + /** + * Return the value of the function + * at the given point. + */ + virtual double operator () (const Point &p) const; + /** + * Return the gradient of the function + * at the given point. + */ + virtual Point gradient (const Point &p) const; +}; + + + + +double RHSPoly<2>::operator () (const Point<2> &p) const { + const double x = p(0), + y = p(1); + const double pi= 3.1415926536; + return 4*pi*pi*(sin(2*pi*x)+sin(2*pi*y)); +}; + + + +double Solution<2>::operator () (const Point<2> &p) const { + const double x = p(0), + y = p(1); + const double pi= 3.1415926536; + return sin(2*pi*x)+sin(2*pi*y); +}; + + +Point<2> Solution<2>::gradient (const Point<2> &p) const { + const double x = p(0), + y = p(1); + const double pi= 3.1415926536; + return Point<2> (2*pi*cos(2*pi*x), + 2*pi*cos(2*pi*y)); +}; + + + + + + +void PoissonEquation<2>::assemble (dFMatrix &cell_matrix, + dVector &rhs, + const FEValues<2> &fe_values, + const Triangulation<2>::cell_iterator &) const { + for (unsigned int point=0; point +void PoissonEquation::assemble (dFMatrix &, + const FEValues &, + const Triangulation::cell_iterator &) const { + Assert (false, ExcPureVirtualFunctionCalled()); +}; + + + +template +void PoissonEquation::assemble (dVector &, + const FEValues &, + const Triangulation::cell_iterator &) const { + Assert (false, ExcPureVirtualFunctionCalled()); +}; + + + + + + + + + +template +PoissonProblem::PoissonProblem (unsigned int order) : + tria(0), dof(0), rhs(0), + boundary_values(0), order(order) {}; + + + + +template +void PoissonProblem::clear () { + if (tria != 0) { + delete tria; + tria = 0; + }; + + if (dof != 0) { + delete dof; + dof = 0; + }; + + if (rhs != 0) + { + delete rhs; + rhs = 0; + }; + + if (boundary_values != 0) + { + delete boundary_values; + boundary_values = 0; + }; + + ProblemBase::clear (); +}; + + + + +template +void PoissonProblem::create_new () { + clear (); + + tria = new Triangulation(); + dof = new DoFHandler (tria); + set_tria_and_dof (tria, dof); +}; + + + + +template +int PoissonProblem::run (const unsigned int level) { + create_new (); + + cout << "Refinement level = " << level + << ", 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 (); + HyperBallBoundary boundary_description; + tria->set_boundary (&boundary_description); + tria->begin_active()->set_refine_flag(); + (++(++(tria->begin_active())))->set_refine_flag(); + tria->execute_refinement (); + tria->refine_global (level); + cout << tria->n_active_cells() << " active cells." << endl; + + rhs = new RHSPoly(); + boundary_values = new Solution (); + + + FiniteElement *fe; + PoissonEquation equation (*rhs); + Quadrature *quadrature; + Quadrature *boundary_quadrature; + switch (order) { + case 0: + fe = new FECrissCross(); + quadrature = new QCrissCross1(); + boundary_quadrature = new QGauss2(); + break; + case 1: + fe = new FELinear(); + quadrature = new QGauss3(); + boundary_quadrature = new QGauss2(); + break; + case 2: + fe = new FEQuadraticSub(); + quadrature = new QGauss4(); + boundary_quadrature = new QGauss3(); + break; + case 3: + fe = new FECubicSub(); + quadrature = new QGauss5(); + boundary_quadrature = new QGauss4(); + break; + case 4: + fe = new FEQuarticSub(); + quadrature = new QGauss6(); + boundary_quadrature = new QGauss5(); + break; + default: + return 100000; + }; + + cout << " Distributing dofs... "; + dof->distribute_dofs (*fe); + cout << dof->n_dofs() << " degrees of freedom." << endl; + n_dofs.push_back (dof->n_dofs()); + + cout << " Assembling matrices..." << endl; + UpdateFlags update_flags = UpdateFlags(update_q_points | update_gradients | + update_JxW_values); + + ProblemBase::FunctionMap dirichlet_bc; + dirichlet_bc[0] = boundary_values; + assemble (equation, *quadrature, *fe, update_flags, dirichlet_bc); + + cout << " Solving..." << endl; + solve (); + + Solution sol; + dVector l1_error_per_cell, l2_error_per_cell, linfty_error_per_cell; + dVector h1_seminorm_error_per_cell, h1_error_per_cell; + + cout << " Calculating L1 error... "; + VectorTools::integrate_difference (*dof_handler, + solution, sol, + l1_error_per_cell, + *quadrature, *fe, L1_norm); + cout << l1_error_per_cell.l1_norm() << endl; + l1_error.push_back (l1_error_per_cell.l1_norm()); + + cout << " Calculating L2 error... "; + VectorTools::integrate_difference (*dof_handler, + solution, sol, + l2_error_per_cell, + *quadrature, *fe, L2_norm); + cout << l2_error_per_cell.l2_norm() << endl; + l2_error.push_back (l2_error_per_cell.l2_norm()); + + cout << " Calculating L-infinity error... "; + VectorTools::integrate_difference (*dof_handler, + solution, sol, + linfty_error_per_cell, + *quadrature, *fe, Linfty_norm); + cout << linfty_error_per_cell.linfty_norm() << endl; + linfty_error.push_back (linfty_error_per_cell.linfty_norm()); + + cout << " Calculating H1-seminorm error... "; + VectorTools::integrate_difference (*dof_handler, + solution, sol, + h1_seminorm_error_per_cell, + *quadrature, *fe, H1_seminorm); + cout << h1_seminorm_error_per_cell.l2_norm() << endl; + h1_seminorm_error.push_back (h1_seminorm_error_per_cell.l2_norm()); + + cout << " Calculating H1 error... "; + VectorTools::integrate_difference (*dof_handler, + solution, sol, + h1_error_per_cell, + *quadrature, *fe, H1_norm); + cout << h1_error_per_cell.l2_norm() << endl; + h1_error.push_back (h1_error_per_cell.l2_norm()); + + if (dof->n_dofs()<=5000) + { + dVector l1_error_per_dof, l2_error_per_dof, linfty_error_per_dof; + dVector h1_seminorm_error_per_dof, h1_error_per_dof; + dof->distribute_cell_to_dof_vector (l1_error_per_cell, l1_error_per_dof); + dof->distribute_cell_to_dof_vector (l2_error_per_cell, l2_error_per_dof); + dof->distribute_cell_to_dof_vector (linfty_error_per_cell, + linfty_error_per_dof); + dof->distribute_cell_to_dof_vector (h1_seminorm_error_per_cell, + 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::project (*dof, constraints, *fe, +// StraightBoundary(), *quadrature, +// sol, projected_solution, false, +// *boundary_quadrature); +// cout << " Calculating L2 error of projected solution... "; +// VectorTools::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; + filename = ('0'+order); + filename += "."; + filename += ('0'+level); + filename += ".ucd"; + cout << " Writing error plots to <" << filename << ">..." << endl; + + DataOut out; + ofstream o(filename.c_str()); + fill_data (out); + 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"); + out.add_data_vector (h1_seminorm_error_per_dof, "H1-seminorm-Error"); + out.add_data_vector (h1_error_per_dof, "H1-Error"); + out.write_ucd (o); + o.close (); + } + else + cout << " Not writing error as grid." << endl; + + cout << endl; + + delete fe; + delete quadrature; + delete boundary_quadrature; + + return dof->n_dofs(); +}; + + +template +void PoissonProblem::print_history (string filename) const { + ofstream out(filename.c_str()); + out << "# n_dofs l1_error l2_error linfty_error h1_seminorm_error h1_error" + << endl; + for (unsigned int i=0; ih/2:" << endl; + cout << " L1 error : " << 1./average_l1 << endl + << " L2 error : " << 1./average_l2 << endl + << " Linfty error : " << 1./average_linfty << endl + << " H1 seminorm error: " << 1./average_h1_semi << endl + << " H1 error : " << 1./average_h1 << endl; + cout << "==========================================================\n"; + cout << "==========================================================\n"; +}; + + + + +int main () { + for (unsigned int order=0; order<5; ++order) + { + PoissonProblem<2> problem (order); + + unsigned int level=0; + unsigned int n_dofs; + do + n_dofs = problem.run (level++); + while (n_dofs<25000); + + string filename; + switch (order) + { + case 0: + filename = "criss_cross"; + break; + case 1: + filename = "linear"; + break; + case 2: + filename = "quadratic"; + break; + case 3: + filename = "cubic"; + break; + case 4: + filename = "quartic"; + break; + }; + filename += ".history"; + + cout << endl << "Printing convergence history to <" + << filename << ">..." << endl; + problem.print_history (filename); + cout << endl << endl << endl; + }; + + return 0; +}; diff --git a/tests/big-tests/multigrid/make_ps b/tests/big-tests/multigrid/make_ps new file mode 100644 index 0000000000..c2df845aec --- /dev/null +++ b/tests/big-tests/multigrid/make_ps @@ -0,0 +1,52 @@ +set term postscript eps +set xlabel "Number of degrees of freedom" +set data style linespoints +set logscale xy + + + +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 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 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 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 output "l2error.eps" +set ylabel "L2-error" + +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 output "h1error.eps" +set ylabel "H1-error" + +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" + + -- 2.39.5