← Back to branch summary

~njansson/dolfin/hpc

« back to all changes in this revision

Viewing changes to demo/ode/homotopy/economy/Economy.h

  • 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:
demo/ode/homotopy/economy/Economy.h
 
1
// Copyright (C) 2005 Anders Logg.
 
2
// Licensed under the GNU LGPL Version 2.1.
 
3
//
 
4
// First added:  2005-04-01
 
5
// Last changed: 2005
 
6
 
 
7
#ifndef __ECONOMY_H
 
8
#define __ECONOMY_H
 
9
 
 
10
#include <dolfin.h>
 
11
 
 
12
using namespace dolfin;
 
13
 
 
14
/// Base class for economies
 
15
 
 
16
class Economy : public Homotopy
 
17
{
 
18
public:
 
19
 
 
20
  Economy(unsigned int m, unsigned int n) :
 
21
    Homotopy(n), m(m), n(n), a(0), w(0),
 
22
    tmp0(0), tmp1(0), tmp2(0), tmp3(0)
 
23
  {
 
24
    a = new real * [m];
 
25
    w = new real * [m];
 
26
    for (unsigned int i = 0; i < m; i++)
 
27
    {
 
28
      a[i] = new real[n];
 
29
      w[i] = new real[n];
 
30
      for (unsigned int j = 0; j < n; j++)
 
31
      {
 
32
        a[i][j] = dolfin::rand();
 
33
        w[i][j] = dolfin::rand();
 
34
      }
 
35
    }
 
36
  }
 
37
 
 
38
  ~Economy()
 
39
  {
 
40
    for (unsigned int i = 0; i < m; i++)
 
41
    {
 
42
      delete [] a[i];
 
43
      delete [] w[i];
 
44
    }
 
45
    delete [] a;
 
46
    delete [] w;
 
47
   
 
48
    if ( tmp0 ) delete [] tmp0;
 
49
    if ( tmp1 ) delete [] tmp1;
 
50
    if ( tmp2 ) delete [] tmp2;
 
51
    if ( tmp3 ) delete [] tmp3;
 
52
  }
 
53
 
 
54
  void disp()
 
55
  {
 
56
    cout << "Trader preferences:" << endl;
 
57
    for (unsigned int i = 0; i < m; i++)
 
58
    {
 
59
      cout << "  trader " << i << ":";
 
60
      for (unsigned int j = 0; j < n; j++)
 
61
        cout << " " << a[i][j];
 
62
      cout << endl;
 
63
    }
 
64
 
 
65
    cout << "Initial endowments:" << endl;
 
66
    for (unsigned int i = 0; i < m; i++)
 
67
    {
 
68
      cout << "  trader " << i << ":";
 
69
      for (unsigned int j = 0; j < n; j++)
 
70
        cout << " " << w[i][j];
 
71
      cout << endl;
 
72
    }
 
73
  }
 
74
 
 
75
  // Number of traders
 
76
  unsigned int m;
 
77
 
 
78
  // Number of goods
 
79
  unsigned int n;
 
80
  
 
81
  // Matrix of traders' preferences
 
82
  real** a;
 
83
  
 
84
  // Matrix of traders' initial endowments
 
85
  real** w;
 
86
 
 
87
protected:
 
88
 
 
89
  // Compute sum of elements
 
90
  complex sum(const complex x[]) const
 
91
  {
 
92
    complex sum = 0.0;
 
93
    for (unsigned int j = 0; j < n; j++)
 
94
      sum += x[j];
 
95
    return sum;
 
96
  }
 
97
 
 
98
  // Compute sum of elements
 
99
  complex bsum(const complex x[], unsigned int b) const
 
100
  {
 
101
    complex sum = 0.0;
 
102
    for (unsigned int j = 0; j < n; j++)
 
103
      sum += std::pow(x[j], b);
 
104
    return sum;
 
105
  }
 
106
 
 
107
  // Compute scalar product x . y
 
108
  complex dot(const real x[], const complex y[]) const
 
109
  {
 
110
    complex sum = 0.0;
 
111
    for (unsigned int j = 0; j < n; j++)
 
112
      sum += x[j] * y[j];
 
113
    return sum;
 
114
  }
 
115
 
 
116
  // Compute special scalar product x . y.^b
 
117
  complex bdot(const real x[], const complex y[], unsigned int b) const
 
118
  {
 
119
    complex sum = 0.0;
 
120
    for (unsigned int j = 0; j < n; j++)
 
121
      sum += x[j] * std::pow(y[j], b);
 
122
    return sum;
 
123
  }
 
124
 
 
125
  // Compute special scalar product x . y.^b
 
126
  complex bdot(const real x[], const complex y[], real b) const
 
127
  {
 
128
    complex sum = 0.0;
 
129
    for (unsigned int j = 0; j < n; j++)
 
130
      sum += x[j] * std::pow(y[j], b);
 
131
    return sum;
 
132
  }
 
133
 
 
134
  // Compute special scalar product x . y.^b
 
135
  complex bdot(const complex x[], const complex y[], unsigned int b) const
 
136
  {
 
137
    complex sum = 0.0;
 
138
    for (unsigned int j = 0; j < n; j++)
 
139
      sum += x[j] * std::pow(y[j], b);
 
140
    return sum;
 
141
  }
 
142
 
 
143
  // Compute special scalar product x . y. z^b
 
144
  complex bdot(const real x[], const complex y[], const complex z[], unsigned int b) const
 
145
  {
 
146
    complex sum = 0.0;
 
147
    for (unsigned int j = 0; j < n; j++)
 
148
      sum += x[j] * y[j] * std::pow(z[j], b);
 
149
    return sum;
 
150
  }
 
151
 
 
152
  // Compute special scalar product x . y. z^b
 
153
  complex bdot(const real x[], const complex y[], const complex z[], real b) const
 
154
  {
 
155
    complex sum = 0.0;
 
156
    for (unsigned int j = 0; j < n; j++)
 
157
      sum += x[j] * y[j] * std::pow(z[j], b);
 
158
    return sum;
 
159
  }
 
160
  
 
161
  // Display values
 
162
  void disp(const real x[], const char* name)
 
163
  {
 
164
    dolfin::cout << name << " = [";
 
165
    for (unsigned int j = 0; j < n; j++)
 
166
      dolfin::cout << x[j] << " ";
 
167
    dolfin::cout << "]" << endl;
 
168
  }
 
169
 
 
170
  // Display values
 
171
  void disp(const complex z[], const char* name)
 
172
  {
 
173
    dolfin::cout << name << " = [";
 
174
    for (unsigned int j = 0; j < n; j++)
 
175
      dolfin::cout << z[j] << " ";
 
176
    dolfin::cout << "]" << endl;
 
177
  }
 
178
 
 
179
  // Initialize temporary storage for scalar products
 
180
  void init(complex** tmp)
 
181
  {
 
182
    *tmp = new complex[m];
 
183
    for (unsigned int i = 0; i < m; i++)
 
184
      (*tmp)[i] = 0.0;
 
185
  }
 
186
  
 
187
  complex* tmp0; // Temporary storage for scalar products
 
188
  complex* tmp1; // Temporary storage for scalar products
 
189
  complex* tmp2; // Temporary storage for scalar products
 
190
  complex* tmp3; // Temporary storage for scalar products
 
191
  
 
192
};
 
193
 
 
194
#endif