~njansson/dolfin/hpc

« back to all changes in this revision

Viewing changes to demo/ode/harmonic/main.cpp

  • Committer: Johannes Ring
  • Date: 2008-03-05 22:43:06 UTC
  • Revision ID: johannr@simula.no-20080305224306-2npsdyhfdpl2esji
The BIG commit!

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
// Copyright (C) 2002 Anders Logg.
 
2
// Licensed under the GNU LGPL Version 2.1.
 
3
//
 
4
// First added:  2002
 
5
// Last changed: 2006-08-22
 
6
//
 
7
// This demo solves the harmonic oscillator on
 
8
// the time interval (0, 4*pi) and computes the
 
9
// error for a set of methods and orders.
 
10
 
 
11
#include <dolfin.h>
 
12
 
 
13
using namespace dolfin;
 
14
 
 
15
class Harmonic : public ODE
 
16
{
 
17
public:
 
18
  
 
19
  Harmonic() : ODE(2, 4.0 * DOLFIN_PI), e(0.0) {}
 
20
 
 
21
  void u0(uBlasVector& u)
 
22
  {
 
23
    u(0) = 0.0;
 
24
    u(1) = 1.0;
 
25
  }
 
26
 
 
27
  void f(const uBlasVector& u, real t, uBlasVector& y)
 
28
  {
 
29
    y(0) = u(1);
 
30
    y(1) = - u(0);
 
31
  }
 
32
 
 
33
  bool update(const uBlasVector& u, real t, bool end)
 
34
  {
 
35
    if ( !end )
 
36
      return true;
 
37
 
 
38
    real e0 = u(0) - 0.0;
 
39
    real e1 = u(1) - 1.0;
 
40
    e = std::max(std::abs(e0), std::abs(e1));
 
41
 
 
42
    return true;
 
43
  }
 
44
  
 
45
  real error()
 
46
  {
 
47
    return e;
 
48
  }
 
49
  
 
50
private:
 
51
 
 
52
  real e;
 
53
 
 
54
};
 
55
 
 
56
int main()
 
57
{
 
58
  dolfin_set("ODE fixed time step", true);
 
59
  dolfin_set("ODE discrete tolerance", 1e-14);
 
60
 
 
61
  for (int q = 1; q <= 5; q++)
 
62
  {
 
63
    dolfin_set("output destination", "silent");
 
64
    dolfin_set("ODE method", "cg");
 
65
    dolfin_set("ODE order", q);
 
66
 
 
67
    Harmonic ode;
 
68
    ode.solve();
 
69
    dolfin_set("output destination", "terminal");
 
70
 
 
71
    message("cG(%d): e = %.3e", q, ode.error());
 
72
  }
 
73
 
 
74
  cout << endl;
 
75
 
 
76
  for (int q = 0; q <= 5; q++)
 
77
  {
 
78
    dolfin_set("output destination", "silent");
 
79
    dolfin_set("ODE method", "dg");
 
80
    dolfin_set("ODE order", q);
 
81
 
 
82
    Harmonic ode;
 
83
    ode.solve();
 
84
    dolfin_set("output destination", "terminal");
 
85
 
 
86
    message("dG(%d): e = %.3e", q, ode.error());
 
87
  }
 
88
 
 
89
  return 0;
 
90
}