← Back to branch summary

~njansson/dolfin/hpc

« back to all changes in this revision

Viewing changes to src/kernel/pde/TimeDependentPDE.cpp

  • Committer: Anders Logg
  • Date: 2007-01-10 09:04:44 UTC
  • mfrom: (1689.1.221 trunk)
  • Revision ID: logg@simula.no-20070110090444-ecyux3n1qnei4i8h
Remove old head

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/kernel/pde/TimeDependentPDE.cpp
 
1
// Copyright (C) 2006 Johan Jansson.
 
2
// Licensed under the GNU GPL Version 2.
 
3
//
 
4
// Modified by Anders Logg 2006.
 
5
//
 
6
// First added:  2006
 
7
// Last changed: 2006-05-04
 
8
 
 
9
//#ifdef HAVE_PETSC_H
 
10
 
 
11
#include <dolfin/dolfin_log.h>
 
12
#include <dolfin/FEM.h>
 
13
#include <dolfin/Matrix.h>
 
14
#include <dolfin/Vector.h>
 
15
#include <dolfin/GMRES.h>
 
16
#include <dolfin/LU.h>
 
17
#include <dolfin/BilinearForm.h>
 
18
#include <dolfin/LinearForm.h>
 
19
#include <dolfin/Mesh.h>
 
20
#include <dolfin/BoundaryCondition.h>
 
21
#include <dolfin/Function.h>
 
22
#include <dolfin/TimeDependentPDE.h>
 
23
#include <dolfin/Parametrized.h>
 
24
#include <dolfin/TimeStepper.h>
 
25
 
 
26
using namespace dolfin;
 
27
 
 
28
TimeDependentPDE::TimeDependentPDE(BilinearForm& a, LinearForm& L, Mesh& mesh, 
 
29
  BoundaryCondition& bc, int N, real T) : GenericPDE(), x(0),
 
30
                                          _a(&a), _Lf(&L),
 
31
                                          _mesh(&mesh), _bc(&bc),
 
32
                                          N(N), t(0), T(T)
 
33
{
 
34
  x = new Vector(N);
 
35
  dotx = new Vector(N);
 
36
 
 
37
  for(unsigned int i = 0; i < this->L().num_functions; i++)
 
38
  {
 
39
    Function* f = this->L().function(i);
 
40
    f->sync(t);
 
41
  }
 
42
}
 
43
//-----------------------------------------------------------------------------
 
44
TimeDependentPDE::~TimeDependentPDE()
 
45
{
 
46
  // Do nothing
 
47
}
 
48
//-----------------------------------------------------------------------------
 
49
dolfin::uint TimeDependentPDE::solve(Function& U)
 
50
{
 
51
  dolfin_assert(_a);
 
52
  dolfin_assert(_Lf);
 
53
  dolfin_assert(_mesh);
 
54
 
 
55
  init(U);
 
56
 
 
57
  // Make sure u is a discrete function associated with the trial space
 
58
//   u.init(*_mesh, _a->trial());
 
59
//   x = &(u.vector());
 
60
//   x = new Vector(N);
 
61
//   dotx = new Vector(N);
 
62
 
 
63
  //int N = FEM::size(*_mesh, _Lf->test());
 
64
 
 
65
  cout << "N: " << N << endl;
 
66
  cout << "x size: " << x->size() << endl;
 
67
 
 
68
 
 
69
  // Initialize ODE (requires x)
 
70
  ode = new TimeDependentODE(*this, N, T);
 
71
  ts = new TimeStepper(*ode);
 
72
 
 
73
  // Write a message
 
74
  dolfin_info("Solving time dependent PDE.");
 
75
 
 
76
//   File  solutionfile("solution.pvd");
 
77
 
 
78
//   solutionfile << U;
 
79
 
 
80
//   int counter = 0;
 
81
 
 
82
  save(U, t);
 
83
 
 
84
  // Start time-stepping
 
85
  while(t < T) {
 
86
//     cout << "t: " << t << endl;
 
87
    
 
88
    preparestep();
 
89
    t = ts->step();
 
90
 
 
91
 
 
92
 
 
93
//     if((counter % 333 * 2) == 0)
 
94
//     {
 
95
//       solutionfile << U;
 
96
//     }
 
97
 
 
98
//     counter++;
 
99
    save(U, t);
 
100
  }
 
101
 
 
102
  return 0;
 
103
}
 
104
//-----------------------------------------------------------------------------
 
105
void TimeDependentPDE::u0(uBlasVector& u)
 
106
{
 
107
  cout << "TimeDependentPDE::u0" << endl;
 
108
 
 
109
  u.copy(*x, 0, 0, u.size());
 
110
}
 
111
//-----------------------------------------------------------------------------
 
112
void TimeDependentPDE::fu(const Vector& x, Vector& dotx, real t)
 
113
{
 
114
}
 
115
//-----------------------------------------------------------------------------
 
116
dolfin::uint TimeDependentPDE::elementdim()
 
117
{
 
118
  dolfin_assert(_a);
 
119
  return _a->trial().elementdim();
 
120
}
 
121
//-----------------------------------------------------------------------------
 
122
BilinearForm& TimeDependentPDE::a()
 
123
{
 
124
  dolfin_assert(_a);
 
125
  return *_a;
 
126
}
 
127
//-----------------------------------------------------------------------------
 
128
LinearForm& TimeDependentPDE::L()
 
129
{
 
130
  dolfin_assert(_Lf);
 
131
  return *_Lf;
 
132
}
 
133
//-----------------------------------------------------------------------------
 
134
Mesh& TimeDependentPDE::mesh()
 
135
{
 
136
  dolfin_assert(_mesh);
 
137
  return *_mesh;
 
138
}
 
139
//-----------------------------------------------------------------------------
 
140
BoundaryCondition& TimeDependentPDE::bc()
 
141
{
 
142
  dolfin_assert(_bc);
 
143
  return *_bc;
 
144
}
 
145
//-----------------------------------------------------------------------------
 
146
void TimeDependentPDE::init(Function& U)
 
147
{
 
148
}
 
149
//-----------------------------------------------------------------------------
 
150
void TimeDependentPDE::save(Function& U, real t)
 
151
{
 
152
}
 
153
void TimeDependentPDE::preparestep()
 
154
{
 
155
}
 
156
//-----------------------------------------------------------------------------
 
157
void TimeDependentPDE::prepareiteration()
 
158
{
 
159
}
 
160
//-----------------------------------------------------------------------------
 
161
TimeDependentODE::TimeDependentODE(TimeDependentPDE& pde, int N, real T) :
 
162
  ODE(N, T), pde(&pde)
 
163
{
 
164
}
 
165
//-----------------------------------------------------------------------------
 
166
void TimeDependentODE::u0(uBlasVector& u)
 
167
{
 
168
  cout << "TimeDependentODE::u0" << endl;
 
169
  pde->u0(u);
 
170
}
 
171
//-----------------------------------------------------------------------------
 
172
void TimeDependentODE::f(const uBlasVector& u, real t, uBlasVector& y)
 
173
{
 
174
  pde->x->copy(u, 0, 0, u.size());
 
175
 
 
176
  pde->prepareiteration();
 
177
 
 
178
  pde->fu(*(pde->x), *(pde->dotx), t);
 
179
 
 
180
  y.copy(*(pde->dotx), 0, 0, u.size());
 
181
}
 
182
//-----------------------------------------------------------------------------
 
183
bool TimeDependentODE::update(const uBlasVector& u, real t, bool end)
 
184
{
 
185
  return true;
 
186
}
 
187
//-----------------------------------------------------------------------------
 
188
 
 
189
//#endif