From: schrage Date: Tue, 23 Mar 1999 14:37:39 +0000 (+0000) Subject: Tutorial/Laplace revisited. Code broken after last changes to deal.II. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=366bb4a23c989357cd730d36f757e52e8cac2bba;p=dealii-svn.git Tutorial/Laplace revisited. Code broken after last changes to deal.II. git-svn-id: https://svn.dealii.org/trunk@1035 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/tutorial/chapter-3.laplace/assemble.html b/deal.II/doc/tutorial/chapter-3.laplace/assemble.html index 31abf12950..b5ca1ef1c4 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/assemble.html +++ b/deal.II/doc/tutorial/chapter-3.laplace/assemble.html @@ -1,13 +1,15 @@ - - DEAL tutorial: the Laplace problem - - - - - + + deal.II tutorial: the Laplace problem + + + + + + +

Assembling the problem

@@ -16,79 +18,75 @@ In order to assemble the matrices we basically need to:

    -
  1. Generate the matrices, i.e. call the DEAL functions that reserve storage space for us. +
  2. Generate the matrices, i.e. call the deal.II + functions that reserve storage space for us.
  3. Calculate the finite element trial functions
  4. Traverse all existing cells and integrate the problem -using the discretized laplace operator +using the discretized Laplace operator
  5. -Traverse all the cell faces and set the appropriate boundary conditions +Insert the local matrices we have used into the global matrix using the +appropriate deal.II functions
  6. -Insert the local matrices we have used into the global matrix using the -appropriate DEAL functions +Traverse all the cell faces and set the appropriate boundary conditions

...and how to do it

-

Function parameters

-
-
-void
-Laplace::assemble_primal(const Function<2>&exact, const Function<2>&)
-{
-
-
-

Generating the matrix structures

-First we generate an n times n square matrix where n is the number -of the degrees of freedom, i.e. the number of points of our discretization. -The parameter max_couplings_between_dofs() returns the maximum -number of couplings between degrees of freedom and allows DEAL -to generate the matrix structure more efficiently, for most of its -elements are zero. +First we generate a structure for the storage of a sparse n times n matrix +where n is the number of the degrees of freedom. +The parameter max_couplings_between_dofs() returns the maximum +number of couplings between degrees of freedom, i.e. the number of elements in the +matrix and allows deal.II to generate the matrix structure more +efficiently, which in effect is done in the next line of code.

-Afterwards the hanging nodes are copied into the matrix, i.e. -the matrix is generated. +Afterwards a constraint matrix for the hanging nodes is created and they +are copied into the matrix structure.

 
-  matrix_structure.reinit(dof_primal.n_dofs(),dof_primal.n_dofs(),
-			  dof_primal.max_couplings_between_dofs());
-  dof_primal.make_sparsity_pattern(matrix_structure);
+  matrix_structure.reinit(dof.n_dofs(),dof.n_dofs(),
+			  dof.max_couplings_between_dofs());
+  dof.make_sparsity_pattern(matrix_structure);
   hanging_nodes.clear();
-  dof_primal.make_constraint_matrix(hanging_nodes);
+  dof.make_constraint_matrix(hanging_nodes);
   hanging_nodes.condense(matrix_structure);
 
 

-The problem is of the form Au=f: +The problem is of the form Au=f, we generate the matrix A with the +structure given by matrix_structure:

 
   A.reinit(matrix_structure);
-  f.reinit(dof_primal.n_dofs());
+  f.reinit(dof.n_dofs());
 
 
-

Calculatinginite element trial functions

+

Calculating finite element trial functions

The two lines below calculate trial functions for the finite elements and -for their faces using Gaussian quadrature. +for their faces using Gaussian quadrature. The first line calculates the trial +function for the finite element associated with the degree of freedom dof +updating the values of the gradients and of the Jacobi determinant multiplied by a +weight function given by the quadrature qc.

 
-  FEValues<2> fevalues(fe_primal, qc_primal, UpdateFlags(update_gradients |
+  FEValues<dim> fevalues(dof->get_fe(), qc, UpdateFlags(update_gradients |
 						   update_JxW_values));
-  FEFaceValues<2> ffvalues(fe_primal, qf_primal,
+  FEFaceValues<dim> ffvalues(dof->get_fe(), qf,
 			   UpdateFlags(update_JxW_values | update_q_points));
 
 
@@ -101,39 +99,55 @@ Integration is done locally. Therefore we need appropriate definitions for
  • an index vector that will allow us to reassemble the global matrix later on
  • -
  • a vector of doubles with the dimension of the total number of degrees of freedom
  • +
  • a vector of doubles with the dimension of the number of degrees of freedom per cell +
  • and a square matrix of doubles with the same dimension
  •  
    -  vector<int> indices(fe_primal.total_dofs);
    -  dVector elvec(fe_primal.total_dofs);
    +  vector<int> indices(fe.total_dofs);
    +  dVector elvec(fe.total_dofs);
       
    -  dFMatrix elmat(fe_primal.total_dofs);
    +  dFMatrix elmat(fe.total_dofs);
     
     

    Next we traverse all the cells and integrate the Laplace problem using the -discretized Laplace operator. qc_primal is a Gaussian quadrature. +discretized Laplace operator. qc is a +Quadrature<dim&rt;.

    +

    +The outer loop traverses all the points of the quadrature qc. +The inner two loops traverse the degrees of freedom of the finite element +fe where du and dv are the gradients +with respect to the quadrature points. fevalues.JxW(k) gives +the Jacobi determinant multiplied by the weight of the quadrature point +k. Taken together the line
    + +elmat(i,j) += fevalues.JxW(k) * du * dv; + +
    +gives the discretized Laplace operator. +

    +
     
    -  for (DoFHandler<2>::active_cell_iterator c = dof_primal.begin_active()
    -					; c != dof_primal.end() ; ++c)
    +  for (DoFHandler<2>::active_cell_iterator c = dof.begin_active()
    +					; c != dof.end() ; ++c)
       {
         fevalues.reinit(c, stb);
         elmat.clear();
         elvec.clear();
    -    c->get_dof_indices(indices);
    +    c->get_dofindices(indices);
         
    -    for (unsigned k=0;k<qc_primal.n_quadrature_points;++k)
    +    for (unsigned k=0;k<qc.n_quadrature_points;++k)
         {
    -      for (unsigned i=0;i<fe_primal.total_dofs;++i)
    +      for (unsigned i=0;i<fe.total_dofs;++i)
           {
     	const Point<2> dv = fevalues.shape_grad(i,k);
     	
    -	for (unsigned j=0;j<fe_primal.total_dofs;++j)
    +	for (unsigned j=0;j<fe.total_dofs;++j)
     	{
     	  const Point<2> du = fevalues.shape_grad(j,k);
     	  
    @@ -146,20 +160,39 @@ discretized Laplace operator. qc_primal is a Gaussian quadrature.
         }
     
     
    +

    +The insertion of the local matrix into the global one +happens in the following piece of code (f is the right hand +vector of the problem, A the problem matrix): +

    +
    +
    +    for (unsigned i=0;i<fe.total_dofs;++i)
    +    {
    +      f(indices[i]) += elvec(i);
    +            
    +      for (unsigned j=0;j<fe.total_dofs;++j)
    +      {
    +	A.add(indices[i], indices[j], elmat(i,j));
    +      }
    +    }
    +
    +

    Setting boundary conditions

    -There are two DEAL functions relevant for us at the moment: +There are two deal.II functions relevant for us at the moment:

     
    -static_void interpolate_boundary_values(...)
    +void VectorTools::interpolate_boundary_values(...)
     
     

    -which does exactly what it says. This function returns a list of pairs -of boundary indicators and the according functions denoting the respective +which does exactly what it says. This function accepts a list of pairs +of boundary indicators and the according functions and returns a list of +pairs of DoF numbers and values denoting the respective Dirichlet boundary values.

    @@ -167,19 +200,19 @@ This output is used by

     
    -static void apply_boundary_values(...)
    +void MatrixTools::apply_boundary_values(...)
     
     

    -that inserts the proper boundary conditions into the equation system: +that inserts the proper boundary conditions into the system of equations:

    
       map<int,double> boundary_values;
       DoFHandler<2>::FunctionMap dirichlet_bc;
    -  BoundaryFct bfkt;
    -  dirichlet_bc[0]=&bfkt;
    -  VectorTools<2>::interpolate_boundary_values(dof_primal,dirichlet_bc,fe_primal,boundary,boundary_values);
    +  BoundaryFct bfct;
    +  dirichlet_bc[0]=&bfct;
    +  VectorTools<2>::interpolate_boundary_values(dof,dirichlet_bc,fe,boundary,boundary_values);
       u.reinit(f);
       MatrixTools<2>::apply_boundary_values(boundary_values,A,u,f);  
     
    @@ -188,11 +221,12 @@ First, we need a few definitions:

    @@ -200,15 +234,15 @@ This may seem a bit confusing. What actually happens is the following:

    1. interpolate_boundary_values takes the boundary functions -bfkt, its relation to boundaries dirichlet_bc and -the triangulation dof_primal, fe_primal and returns a +bfct, its relation to boundaries dirichlet_bc and +the triangulation dof, fe and returns a mapping boundary_values that maps values instead of functions to our boundaries. The function looks at all the boundaries. All we ever need to do is specify the initial triangulation.
    2. apply_boundary_values subsequently takes that mapping and -our equation system Au=f and inserts the boundary values into -the equation system which can then be solved. +our system of equations Au=f and inserts the boundary values into +the system of equations which can then be solved.
    @@ -221,7 +255,7 @@ the equation system which can then be solved.
    Jan Schrage

    -Last modified: Fri Feb 12, 1999 +Last modified: Tue Mar 9, 1999

    diff --git a/deal.II/doc/tutorial/chapter-3.laplace/code/func.cc b/deal.II/doc/tutorial/chapter-3.laplace/code/func.cc index a5c558ea79..8a739c6b2e 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/code/func.cc +++ b/deal.II/doc/tutorial/chapter-3.laplace/code/func.cc @@ -1,6 +1,5 @@ // $Id$ -// JS. const char* funcversion = "Functions: $Revision$"; #include "functions.h" diff --git a/deal.II/doc/tutorial/chapter-3.laplace/code/functions.h b/deal.II/doc/tutorial/chapter-3.laplace/code/functions.h index d8b592d352..9b888fd69e 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/code/functions.h +++ b/deal.II/doc/tutorial/chapter-3.laplace/code/functions.h @@ -1,6 +1,5 @@ // $Id$ -// JS.Wird das File ueberhaupt gebraucht ? #include diff --git a/deal.II/doc/tutorial/chapter-3.laplace/code/laplace.cc b/deal.II/doc/tutorial/chapter-3.laplace/code/laplace.cc index 0b4d0d5f56..da588e58ac 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/code/laplace.cc +++ b/deal.II/doc/tutorial/chapter-3.laplace/code/laplace.cc @@ -20,120 +20,107 @@ const char* laplaceversion = "Laplace: $Revision$"; #include #include #include - -#include +#include +#include #include #include #include #define PRIMEL FELinear<2> -#define DUEL FEQuadraticSub<2> // Finite Elements -static PRIMEL fe_primal; -// JS.static DUEL fe_dual; +static PRIMEL fe; // Quadrature formulae -static QGauss2<2> qc_primal; -static QGauss2<1> qf_primal; -// static QGauss3<2> qc_dual; -// static QGauss3<1> qf_dual; +static QGauss2<2> qc; +static QGauss2<1> qf; static QGauss5<2> qc_integrate; static QGauss5<1> qf_integrate; StraightBoundary<2> stb; -// JS.ist im Moment noch PureTransportSolution... Laplace::Laplace() - : dof_primal(&tr) + : dof(&tr) { - // JS.Triangulation generieren. Zellränder werden numeriert. + // Generate the initial triangulation tr.create_hypercube(-1.,1.); - // JS.Freiheitsgrade verteilen. D.h. Zellen numerieren. - dof_primal.distribute_dofs(fe_primal); + // Distribute the degrees of freedom + dof.distribute_dofs(fe); } Laplace::~Laplace() {} -// JS.Gitter verfeinern. +// Remesh the grid void Laplace::remesh(unsigned int steps) { if (tr.n_levels() <= 1) { - tr.refine_global(1); //JS.Lokal ist execute_coarsening_etc... + tr.refine_global(1); // refine globally } if (steps) tr.refine_global(steps); else - tr.execute_coarsening_and_refinement(); - - // JS. Freiheitsgrade neu verteilen. - dof_primal.distribute_dofs(fe_primal); - // JS. und dem Problem angemessener nochmal numerieren. - dof_primal.renumber_dofs(Cuthill_McKee); - deallog << "Cells " << tr.n_active_cells() - << " PrimalDoFs " << dof_primal.n_dofs() - << endl; + tr.execute_coarsening_and_refinement(); // refine locally + + // redistribute the degrees of freedom... + dof.distribute_dofs(fe); + // ...and renumber them so they can be used more efficiently + dof.renumber_dofs(Cuthill_McKee); } // JS.Primales Problem zusammenstellen. void -Laplace::assemble_primal() +Laplace::assemble() { - deallog << "Assembling primal problem" << endl; - // JS. Platz für neue Matrix mit (?) (2x = quadratisch) Anzahl der Zellen, - // Anzahl der Kopplungen (Matrix dünn besetzt, für effizientes Speichern) - matrix_structure.reinit(dof_primal.n_dofs(),dof_primal.n_dofs(), - dof_primal.max_couplings_between_dofs()); - // JS.Hängende Noden in die Matrix einbauen; d.h. Matrix generieren. - dof_primal.make_sparsity_pattern(matrix_structure); + // Initialize the problem matrix, i.e. reserve storage space + matrix_structure.reinit(dof.n_dofs(),dof.n_dofs(), + dof.max_couplings_between_dofs()); + // Generate the matrix structure + dof.make_sparsity_pattern(matrix_structure); hanging_nodes.clear(); - dof_primal.make_constraint_matrix(hanging_nodes); + dof.make_constraint_matrix(hanging_nodes); hanging_nodes.condense(matrix_structure); - //JS.Problem der Form Au=f. + // The problem has the form Au=f. A.reinit(matrix_structure); - f.reinit(dof_primal.n_dofs()); + f.reinit(dof.n_dofs()); - // JS.Ansatzfunktionen auf Zellrändern im Voraus berechnen aus - // Effizienzgründen. - FEValues<2> fevalues(fe_primal, qc_primal, UpdateFlags(update_gradients | + // Calculate the trial functions on the cell faces. + FEValues<2> fevalues(fe, qc, UpdateFlags(update_gradients | update_JxW_values)); - FEFaceValues<2> ffvalues(fe_primal, qf_primal, + FEFaceValues<2> ffvalues(fe, qf, UpdateFlags(update_JxW_values | update_q_points)); - //JS.Ab hier lokales Problem el... = Finites Element... - //JS. Index für eine Zelle, für späteren Einbau in globale Matrix. - vector indices(fe_primal.total_dofs); - dVector elvec(fe_primal.total_dofs); + + // Integrate the problem locally... + vector indices(fe.total_dofs); + Vector elvec(fe.total_dofs); - dFMatrix elmat(fe_primal.total_dofs); + FullMatrix elmat(fe.total_dofs); - // JS.Einmal alle Zellen durchlaufen - for (DoFHandler<2>::active_cell_iterator c = dof_primal.begin_active() - ; c != dof_primal.end() ; ++c) + for (DoFHandler<2>::active_cell_iterator c = dof.begin_active() + ; c != dof.end() ; ++c) { fevalues.reinit(c, stb); elmat.clear(); elvec.clear(); c->get_dof_indices(indices); - // JS.Integration des Problems. Diese Schleifenfolge für Effizienz. - for (unsigned k=0;k dv = fevalues.shape_grad(i,k); - for (unsigned j=0;j du = fevalues.shape_grad(j,k); @@ -144,74 +131,56 @@ Laplace::assemble_primal() } } } - // JS.Lokale Matrix in globale einbauen. - for (unsigned i=0;i boundary_values; DoFHandler<2>::FunctionMap dirichlet_bc; BoundaryFct bfkt; dirichlet_bc[0]=&bfkt; - VectorTools<2>::interpolate_boundary_values(dof_primal,dirichlet_bc, - fe_primal,boundary, + VectorTools<2>::interpolate_boundary_values(dof,dirichlet_bc, + fe,boundary, boundary_values); u.reinit(f); MatrixTools<2>::apply_boundary_values(boundary_values,A,u,f); - - - // JS.Hängende Noden einbauen. - // hanging_nodes.condense(A); - //hanging_nodes.condense(f); } -// JS. Primales Problem lösen. +// Solve the primal problem. void -Laplace::solve_primal() +Laplace::solve() { - deallog.push("Solve"); - - // JS.Empfindlichkeit des Lösers einstellen. + // Solver control: max 1000 iterations, threshold 1e-10 SolverControl control(1000, 1.e-10); - // JS. Löser definieren. modifiziertes cg-Verfahren, - SolverCG solver(control, mem); - - // JS.??? - // u.reinit(f); + // Define the solver. + SolverCG> solver(control, mem); - // JS.lösen. solver.solve(A,u,f); - // JS.??? hanging_nodes.distribute(u); - - deallog.pop(); } -// JS. Datenausgabe im Gnuplot-Format. +// Data output for gnuplot. void Laplace::write_data(const char* name) { - deallog << "Writing gnuplot" << endl; DataOut<2> out; char fname[100]; - sprintf(fname,"P_%s",name); - { ofstream gnuplot(fname); out.clear_data_vectors(); - out.attach_dof_handler(dof_primal); + out.attach_dof_handler(dof); out.add_data_vector(u,"solution","kg"); out.write_gnuplot (gnuplot, 1); diff --git a/deal.II/doc/tutorial/chapter-3.laplace/code/laplace.h b/deal.II/doc/tutorial/chapter-3.laplace/code/laplace.h index 00cc5a017a..552ea2a101 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/code/laplace.h +++ b/deal.II/doc/tutorial/chapter-3.laplace/code/laplace.h @@ -7,16 +7,16 @@ #include #include #include -#include -#include +#include +#include -class AdvMatrix : - public dSMatrix +class AdvMatrix : + public SparseMatrix { public: - void precondition(dVector& dst, const dVector& src) const + void precondition(Vector& dst, const Vector& src) const { - dSMatrix::precondition_SSOR(dst, src); + SparseMatrix::precondition_SSOR(dst, src); } }; @@ -26,16 +26,16 @@ class Laplace protected: Point<2> direction; Triangulation<2> tr; - DoFHandler<2> dof_primal; + DoFHandler<2> dof; - dSMatrixStruct matrix_structure; + SparseMatrixStruct matrix_structure; AdvMatrix A; - dVector u; - dVector z; - dVector f; + Vector u; + Vector z; + Vector f; - PrimitiveVectorMemory mem; + PrimitiveVectorMemory mem; ConstraintMatrix hanging_nodes; @@ -46,8 +46,8 @@ public: ~Laplace(); void remesh(unsigned int global_refine = 0); - void assemble_primal(); - void solve_primal(); + void assemble(); + void solve(); void write_data(const char* name); }; diff --git a/deal.II/doc/tutorial/chapter-3.laplace/code/main.cc b/deal.II/doc/tutorial/chapter-3.laplace/code/main.cc index 73a8900e24..76eb23f5eb 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/code/main.cc +++ b/deal.II/doc/tutorial/chapter-3.laplace/code/main.cc @@ -7,43 +7,23 @@ #include #include -#include -#include - -ZeroFunction<2> zero; - char fname[30]; main(int argc, char** argv) { - // JS.Logfile erzeugen, als Stream konzipiert - - ofstream logfile("T"); - deallog.attach(logfile); - if (argc==1) cerr << "Usage: " << argv[0] << "firstgrid\nUsing 3" << endl; int firstgrid = 3; - - if (argc>=2) firstgrid = atoi(argv[1]); - deallog << "Firstgrid " << firstgrid << endl; + if (argc>=2) firstgrid = atoi(argv[1]); - // JS.Benötigte Funktionen zur Lösung des Problems Laplace lap; - // JS.Logstream ist ein stack, ab hier wird "Adaptive:" vor jede Zeile - // gestellt - deallog.push("Adaptive"); - deallog.depth_console(2); - for (unsigned step = 0; step < 3 ; ++step) { - deallog << "Step " << step << endl; - // JS.Beim ersten Mal ein verfeinertes Grid erzeugen, firstgrid=Anzahl - // der Verfeinerungen + // Generate a refined grid. if (!step) lap.remesh(firstgrid); else @@ -51,20 +31,12 @@ main(int argc, char** argv) lap.remesh(1); } - deallog.push("Primal"); - - // JS.exakte Lösung mit 0 auf der rechten Seite - lap.assemble_primal(); - lap.solve_primal(); - - // JS.ab hier kein "Adaptive:" mehr - deallog.pop(); - sprintf(fname,"T%02d",step); - // JS.Daten zurückschreiben. + // Assemble and solve the problem. + lap.assemble(); + lap.solve(); + + // Data output. lap.write_data(fname); } - // JS.Logfile sauber abschließen und Schluß - deallog.pop(); - deallog.detach(); } diff --git a/deal.II/doc/tutorial/chapter-3.laplace/index.html b/deal.II/doc/tutorial/chapter-3.laplace/index.html index da40543860..39d39b01b6 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/index.html +++ b/deal.II/doc/tutorial/chapter-3.laplace/index.html @@ -3,17 +3,16 @@ The Laplace Problem - -i - + + +

    The Laplace Problem

    We will start with a simple example for a differential equation to be solved -with DEAL: the Laplace problem. We will try to solve the +with deal.II: the Laplace problem. We will try to solve the Laplace equation on a square where one of the four boundaries has a constant and the other three zero potential.

    @@ -31,7 +30,7 @@ laplace problem assemble matrices describing the boundary conditions

  • -let DEAL solve the problem +let deal.II solve the problem
  • take care of the data output @@ -50,7 +49,7 @@ where the overall design of the program and its classes are discussed

  • The main program

    -where we take a look at how DEAL is used

    +where we take a look at how deal.II is used

  • The class Laplace @@ -58,11 +57,9 @@ where we take a look at how DEAL is used

    where some of the details of the problem generation and solution are discussed

  • -
  • -Generating a -triangulation

    -

    -where a triangulation is generated and degrees of freedom are discussed

    +
  • +Generating a triangulation +

    where a triangulation is generated and degrees of freedom are discussed

  • Assembling the problem matrix @@ -73,7 +70,15 @@ and the boundary conditions are set
  • Solving the problem +

    where the problem is solved +

    +
  • +
  • +The source code +

    +contains link to all the source files. +


  • @@ -85,7 +90,7 @@ where the problem is solved
    Jan Schrage

    -Last modified: Mon 15 Feb 1999 +Last modified: Tue 9 Mar 1999

    diff --git a/deal.II/doc/tutorial/chapter-3.laplace/laplace.html b/deal.II/doc/tutorial/chapter-3.laplace/laplace.html index bbc49f1be2..cce0c9f999 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/laplace.html +++ b/deal.II/doc/tutorial/chapter-3.laplace/laplace.html @@ -22,10 +22,9 @@ Let's have a look at the class definition: class Laplace { - Function<2>& exact; -protected: + Triangulation<2> tr; - DoFHandler<2> dof_primal; + DoFHandler<2> dof; dSMatrixStruct matrix_structure; LapMatrix A; @@ -35,11 +34,11 @@ protected:

    -These few lines define several important elements: The right hand side of the equation -exact, the triangulation tr, i.e. the grid, and a handler for the degrees -of freedom for the finite elements dof_primal, all for the two-dimensional case. -In addition three matrices are defined (the matrix A defining our problem). Note that -in order to solve any problem at all with DEAL the definitions above are paramount. +These few lines define several important elements: The triangulation tr, i.e. the grid, and a handler for the degrees +of freedom for the finite elements dof, all for the two-dimensional case. +In addition one matrix structure and two matrices are defined (the matrix A defining our problem). Further explanation can be found in the +chapter Assembling the problem. Note that +in order to solve any problem at all with deal.II the definitions above are paramount.

    @@ -49,31 +48,19 @@ The constructor has the task of generating a triangulation, too.

     
     public:
    -  Laplace(Function<2>& solution);
    +  Laplace();
       ~Laplace();
     
     

    -The next few functions refine the grid - non-adaptively - assemble the primal problem and +The next few functions refine the grid - non-adaptively - assemble the problem and call the appropriate solver.

     
       void remesh(unsigned int global_refine = 0);
    -  void assemble_primal(const Function<2>& boundary, const Function<2>& rhs);
    -  void solve_primal();
    -
    -
    - + void assemble(); + void solve(); };
    diff --git a/deal.II/doc/tutorial/chapter-3.laplace/main.html b/deal.II/doc/tutorial/chapter-3.laplace/main.html index 8093c61e65..a53c8dca34 100644 --- a/deal.II/doc/tutorial/chapter-3.laplace/main.html +++ b/deal.II/doc/tutorial/chapter-3.laplace/main.html @@ -12,7 +12,7 @@

    The Laplace Problem

    -

    Main Program: main.cc

    +

    Main Program: main.cc

    The main program has several functions:

    @@ -46,8 +46,8 @@ the function definition for our solution function are included:
     
    -#include "laplace.h"
    -#include "functions.h"
    +#include "laplace.h"
    +#include "functions.h"
     
     
    @@ -78,7 +78,10 @@ Laplace lap;

    In addition we need to do some refinement (the command line argument was -previously stored in the variable firstgrid. +previously stored in the variable firstgrid. During the first +execution +of the loop global refinement is done firstgrid times, every +following time it is done once more.

    
     for (unsigned step = 0; step < 3 ; ++step)   
    @@ -94,13 +97,11 @@ for (unsigned step = 0; step < 3 ; ++step)
     
     

    Problem assemblage and solution

    -Our class assembles and solves the primal problem; the solution is exact (as defined above) -and the right hand side of the equation is zero (as defined above). If the right -hand side were not zero we would solve the Poisson equation instead. +Our class assembles and solves the primal problem.

    
    -  lap.assemble_primal(boundary, zero); 
    -  lap.solve_primal();          
    +  lap.assemble(); 
    +  lap.solve();          
     

    Data output

    @@ -108,7 +109,6 @@ hand side were not zero we would solve the Poisson equation instead. Finally the solution is written to a file.

     
    -  sprintf(fname,"T%02d",step);    
       lap.write_data(fname);   
      }   
     }
    diff --git a/deal.II/doc/tutorial/chapter-3.laplace/solution.html b/deal.II/doc/tutorial/chapter-3.laplace/solution.html
    index 5437cd91a8..d728f9fad5 100644
    --- a/deal.II/doc/tutorial/chapter-3.laplace/solution.html
    +++ b/deal.II/doc/tutorial/chapter-3.laplace/solution.html
    @@ -11,6 +11,28 @@
     
     

    Solving the problem

    +

    +Now that the problem matrix is assembled we need to define an appropriate +solver. In this case we use a CG-solver with a maximum of 1000 +iterations and a threshold of 1e-10. +

    + +
    
    +void
    +Laplace::solve_primal()
    +{
    +
    +  SolverControl control(1000, 1.e-10);
    +
    +  SolverCG solver(control, mem);
    +  
    +  solver.solve(A,u,f);
    +  
    +  hanging_nodes.distribute(u);
    +
    +}
    +
    +

    Back to the tutorial index diff --git a/deal.II/doc/tutorial/chapter-3.laplace/source.html b/deal.II/doc/tutorial/chapter-3.laplace/source.html new file mode 100644 index 0000000000..c55dbbbe81 --- /dev/null +++ b/deal.II/doc/tutorial/chapter-3.laplace/source.html @@ -0,0 +1,46 @@ + + + +The Laplace Problem + +i + + + +

    The Source Code for the Solution of the Laplace Problem

    + +You can have a look at the complete source code for the solution of the +Laplace problem by following the links below: +
      +
    • main.cc containing the main program which + sets some parameters and starts the solution process. +
    • +
    • laplace.h and + laplace.cc which define the class + Laplace that assembles the laplace problem. +
    • +
    • functions.h and + func.cc which define the function that sets + the boundary conditions. +
    • +
    • Makefile. Last but not least. You will + probably need to edit the variable root defined at the + top of the Makefile. This defines the location of your + DEAL directory. +
    • +
    + +
    +

    +Back to the tutorial index +

    +
    +
    +Jan Schrage
    +

    +Last modified: Mon 8 Mar 1999 +

    + +