← Back to branch summary

~ubuntu-branches/debian/squeeze/gmsh/squeeze

« back to all changes in this revision

Viewing changes to Numeric/FunctionSpace.h

  • Committer: Bazaar Package Importer
  • Author(s): Christophe Prud'homme, Christophe Prud'homme
  • Date: 2009-09-27 17:36:40 UTC
  • mfrom: (1.2.9 upstream) (8.1.2 squeeze)
  • Revision ID: james.westby@ubuntu.com-20090927173640-oxyhzt0eadjfrlwz
[Christophe Prud'homme]
* New upstream release
  + solver code refactoring
  + better IDE integration.

Show diffs side-by-side

added added

removed removed

Lines of Context:
Numeric/FunctionSpace.h
1
 
// Gmsh - Copyright (C) 1997-2009 C. Geuzaine, J.-F. Remacle
2
 
//
3
 
// See the LICENSE.txt file for license information. Please report all
4
 
// bugs and problems to <gmsh@geuz.org>.
5
 
 
6
 
#ifndef _FUNCTION_SPACE_H_
7
 
#define _FUNCTION_SPACE_H_
8
 
 
9
 
#include <math.h>
10
 
#include <map>
11
 
#include "GmshMatrix.h"
12
 
 
13
 
struct gmshFunctionSpace 
14
 
{
15
 
  gmshMatrix<double> points;
16
 
  gmshMatrix<double> monomials;
17
 
  gmshMatrix<double> coefficients;
18
 
  inline void evaluateMonomials(double u, double v, double w, double p[]) const 
19
 
  {
20
 
    for (int j = 0; j < monomials.size1(); j++) {
21
 
      p[j] = pow(u, (int)monomials(j, 0));
22
 
      if (monomials.size2() > 1) p[j] *= pow(v, (int)monomials(j, 1));
23
 
      if (monomials.size2() > 2) p[j] *= pow(w, (int)monomials(j, 2));
24
 
    }
25
 
  }
26
 
  inline void f(double u, double v, double w, double *sf) const
27
 
  {
28
 
    double p[256];
29
 
    evaluateMonomials(u, v, w, p);
30
 
    for (int i = 0; i < coefficients.size1(); i++) {
31
 
      sf[i] = 0;
32
 
      for (int j = 0; j < coefficients.size2(); j++) {
33
 
        sf[i] += coefficients(i, j) * p[j];
34
 
      }
35
 
    }
36
 
  }
37
 
  inline void df(double u, double v, double w, double grads[][3]) const
38
 
  {
39
 
    switch (monomials.size2()) {
40
 
    case 1:
41
 
      for (int i = 0; i < coefficients.size1(); i++){
42
 
        grads[i][0] = 0;
43
 
        grads[i][1] = 0;
44
 
        grads[i][2] = 0;
45
 
        for(int j = 0; j < coefficients.size2(); j++){
46
 
          if ((monomials)(j, 0) > 0)
47
 
            grads[i][0] += (coefficients)(i, j) * 
48
 
              pow(u, (monomials)(j, 0) - 1) * (monomials)(j, 0);
49
 
        }
50
 
      }
51
 
      break;
52
 
    case 2:
53
 
      for (int i = 0; i < coefficients.size1(); i++){
54
 
        grads[i][0] = 0;
55
 
        grads[i][1] = 0;
56
 
        grads[i][2] = 0;
57
 
        for(int j = 0; j < coefficients.size2(); j++){
58
 
          if ((monomials)(j, 0) > 0)
59
 
            grads[i][0] += (coefficients)(i, j) *
60
 
              pow(u, (monomials)(j, 0) - 1) * (monomials)(j, 0) *
61
 
              pow(v, (monomials)(j, 1));
62
 
          if ((monomials)(j, 1) > 0)
63
 
            grads[i][1] += (coefficients)(i, j) *
64
 
              pow(u, (monomials)(j, 0)) *
65
 
              pow(v, (monomials)(j, 1) - 1) * (monomials)(j, 1);
66
 
        }
67
 
      }
68
 
      break;
69
 
    case 3:
70
 
      for (int i = 0; i < coefficients.size1(); i++){
71
 
        grads[i][0] = 0;
72
 
        grads[i][1] = 0;
73
 
        grads[i][2] = 0;
74
 
        for(int j = 0; j < coefficients.size2(); j++){
75
 
          if ((monomials)(j, 0) > 0)
76
 
            grads[i][0] += (coefficients)(i, j) *
77
 
              pow(u, (monomials)(j, 0) - 1) * (monomials)(j, 0) *
78
 
              pow(v, (monomials)(j, 1)) *
79
 
              pow(w, (monomials)(j, 2));
80
 
          if ((monomials)(j, 1) > 0)
81
 
            grads[i][1] += (coefficients)(i, j) *
82
 
              pow(u, (monomials)(j, 0)) *
83
 
              pow(v, (monomials)(j, 1) - 1) * (monomials)(j, 1) *
84
 
              pow(w, (monomials)(j, 2));
85
 
          if ((monomials)(j, 2) > 0)
86
 
            grads[i][2] += (coefficients)(i, j) *
87
 
              pow(u, (monomials)(j, 0)) *
88
 
              pow(v, (monomials)(j, 1)) *
89
 
              pow(w, (monomials)(j, 2) - 1) * (monomials)(j, 2);
90
 
        }
91
 
      }
92
 
      break;
93
 
    }
94
 
  }
95
 
};
96
 
 
97
 
class gmshFunctionSpaces 
98
 
{
99
 
 private:
100
 
  static std::map<int, gmshFunctionSpace> fs;
101
 
  static std::map<std::pair<int, int>, gmshMatrix<double> > injector;
102
 
 public :
103
 
  static const gmshFunctionSpace &find(int);
104
 
  static const gmshMatrix<double> &findInjector(int, int);
105
 
};
106
 
 
107
 
#endif